Nature | www.nature.com | 1 Article Polarized thermal emission from dust in a galaxy at redshift 2.6 J. E. Geach1 ✉, E. Lopez-Rodriguez2, M. J. Doherty1, Jianhang Chen3, R. J. Ivison3,4,5,6, G. J. Bendo7, S. Dye8 & K. E. K. Coppin1 Magnetic fields are fundamental to the evolution of galaxies, playing a key role in the astrophysics of the interstellar medium and star formation. Large-scale ordered magnetic fields have been mapped in the Milky Way and nearby galaxies1,2, but it is not known how early in the Universe such structures formed3. Here we report the detection of linearly polarized thermal emission from dust grains in a strongly lensed, intrinsically luminous galaxy that is forming stars at a rate more than 1,000 times that of the Milky Way at redshift 2.6, within 2.5 Gyr of the Big Bang4,5. The polarized emission arises from the alignment of dust grains with the local magnetic field6,7. The median polarization fraction is of the order of 1%, similar to nearby spiral galaxies8. Our observations support the presence of a 5-kiloparsec-scale ordered magnetic field with a strength of around 500 μG or lower, oriented parallel to the molecular gas disk. This confirms that such structures can be rapidly formed in galaxies, early in cosmic history. We observed the lensed galaxy 9io9 (ref. 4) with the Atacama Large Millimeter/Submillimeter Array (ALMA) at a representative frequency of 242 GHz (equivalent to a wavelength of roughly 350 μm in the rest- frame of the galaxy) to record the dust continuum emission averaged over a total bandwidth of 7.5 GHz. The set of XX, YY, XY and YX linear polarization parameters recorded in full polarization mode allow meas- urement of the Stokes parameters Q and U, yielding the total linearly polarized intensity, Q UPI = +2 2 and position angle of polarized emission χ U Q= 0.5 arctan( / ). The root mean squared sensitivity of the observations is σI = 47 and σQ ≅ σU = 9 μJy beam−1. In Fig. 1 we present image plane maps of the total intensity I, Stokes Q and U, and polarized intensity (PI). The polarization angle χ is rotated by 90° to show the plane-of-the-sky magnetic (B) field orientation (χB). We measure an image plane integrated flux density of I = 62 mJy, integrated polariza- tion fraction of P = 0.6 ± 0.1%, where P = PI/I and B field orientation of χB = (−0.7 ± 1.4)°. The mean of the distribution of polarization fractions and B field orientations is ⟨P⟩ = 0.6 ± 0.3% and ⟨χB⟩ = (0.8 ± 18.3)°, respec- tively. Note that the uncertainties are the dispersion of the distribution of individual measurements within the galaxy, not the accuracy in the polarization measurement (Methods). Using the lens model derived from previous high-resolution milli- metre continuum emission and optical Hubble Space Telescope imag- ing, the source plane CO(4–3) emission, tracing the cold molecular gas reservoir, has been shown to be well modelled by a rotating disk of maximum radius 2.6 kpc, inclined by roughly 50° to the line of sight, with a position angle on the sky of roughly 5° east of north4,5. With this model as a constraint, we explore what source plane B field configura- tions are consistent with the image plane polarization observations. The most likely source plane configuration is a large-scale ordered B field oriented χ = 5B −10 +5° east of north with an extent matching that of the CO emission (Fig. 2). This result implies the presence of a 5-kpc-scale galactic ordered magnetic field oriented parallel to the molecular gas-rich disk. Angular variations of χB across the galaxy present in the image plane maps, corresponding to scales of 600 pc in the source plane, can be explained by the low signal to noise ratio and beam effects (Methods). This result indicates that the introduction of a random B field component with an angular variation of ±5° in addition to the large-scale ordered B field is also consistent with the observations. A present, we lack the sensitivity and resolution to map the configura- tion of the B field at scales of roughly 100 pc in which structure related to turbulence can start to be resolved. The observed B field configura- tion parallel to the disk is consistent with the galactic B fields measured in local spiral galaxies observed at far-infrared and radio wavelengths1,2. Note that our far-infrared polarimetric observations trace a density- weighted average B field in the cold and dense interstellar medium (ISM), rather than a volume average B field in the warm and diffuse ISM by radio polarimetric observations. The mean and integrated polarization fractions of 9io9 are consistent with the P ≅ 0.8% level measured in nearby spiral and starburst galaxies at wavelengths of 53–214 μm (ref. 2). The observations presented here are sensitive to polarized emission beyond this range, pushing into the Rayleigh–Jeans tail of the thermal emission spectrum at a rest-frame wavelength of λrest = 350 μm. In recent models of diffuse interstellar dust, the polarization fraction, P, is independent of wavelength across 200–2,000 μm, consistent with observations of Galactic dust emission. Observations of local starburst galaxies show that P only varies by 0.4% over the 50–150 μm range, with an increase of up to roughly 1% towards 214 μm (ref. 2). We therefore conclude that 9io9 has a polarization https://doi.org/10.1038/s41586-023-06346-4 Received: 17 February 2023 Accepted: 20 June 2023 Published online: xx xx xxxx Open access Check for updates 1Centre for Astrophysics Research, School of Physics, Engineering and Computer Science, University of Hertfordshire, Hatfield, UK. 2Kavli Institute for Particle Astrophysics and Cosmology, Stanford University, Stanford, CA, USA. 3European Southern Observatory, Garching, Germany. 4Department of Physics and Astronomy, Macquarie University, Sydney, New South Wales, Australia. 5School of Cosmic Physics, Dublin Institute for Advanced Studies, Dublin, Ireland. 6Institute for Astronomy, Royal Observatory, University of Edinburgh, Edinburgh, UK. 7UK ALMA Regional Centre Node, Jodrell Bank Centre for Astrophysics, Department of Physics and Astronomy, The University of Manchester, Manchester, UK. 8School of Physics and Astronomy, University of Nottingham, Nottingham, UK. ✉e-mail: j.geach@herts.ac.uk 2 | Nature | www.nature.com Article level similar to local star-forming disks and starburst galaxies, with a key difference being the order-of-magnitude difference in gas mass and star-formation rate, with the disk of 9io9 being close to molecular gas dominated, contrasted with the fgas ≅ 10% gas fractions of local star-forming disks9. The large-scale ordered magnetic fields that exist in massive disk galaxies in the local Universe is thought to arise through the amplifica- tion of seed fields, and this has been predicted to occur on relatively short cosmological timescales, of order 1 Gyr (refs. 10–12). Weak seed fields (as low as B ≅ 10−20 G) could be formed in protogalaxies either through trapping of a cosmological field, possibly primordial in nature, or through the battery effect following the onset of star formation13–16. Although turbulent gas motions in disks can reduce net polarization if they impart a strong turbulent component to the B field17, recent theo- retical models of the formation of galactic-scale magnetic fields invoke turbulence in the ISM as the origin of a ‘small-scale’ dynamo that can rapidly amplify the weak seed fields to μG levels12,18,19. This small-scale dynamo is mainly driven by supernova explosions with coherence lengths of order 50–100 pc, but turbulence can be injected into the ISM on multiple scales through disk instabilities and feedback effects, includ- ing stellar winds and outflows driven by radiation pressure, supernova explosions and large-scale outflows from an active galactic nucleus. The average turbulent velocity component of the disk of 9io9, deter- mined from kinematic modelling of the CO emission, is σv ≅ 70 km s−1 0° 16′ 02″ 00″ 15′ 58″ 56″ D ec . ( IC R S ) P = 1% 0 2 4 6 8 10 I (mJy per beam) −0.04 −0.02 0 0.02 0.04 −0.04 −0.02 0 0.02 0.04 P = 1% 0 0.02 0.04 0.06 0.08 0.10 02 h 0 9 m in 41 .5 s 41 .4 s 41 .3 s 41 .2 s 41 .1 s 0° 16′ 02″ 00″ 15′ 58″ 56″ RA (J2000) D ec . ( J2 00 0) P = 1% RA (J2000) RA (J2000) RA (J2000) P = 1% a b c d e f g h Stokes I Stokes Q Stokes U PI Q (mJy per beam) 0° 16′ 02″ 00″ 15′ 58″ 56″ 0° 16′ 02″ 00″ 15′ 58″ 56″ 0° 16′ 02″ 00″ 15′ 58″ 56″ 0° 16′ 02″ 00″ 15′ 58″ 56″ U (mJy per beam) 0° 16′ 02″ 00″ 15′ 58″ 56″ 0° 16′ 02″ 00″ 15′ 58″ 56″ PI (mJy per beam) 02 h 0 9 m in 41 .5 s 41 .4 s 41 .3 s 41 .2 s 41 .1 s 02 h 0 9 m in 41 .5 s 41 .4 s 41 .3 s 41 .2 s 41 .1 s 02 h 0 9 m in 41 .5 s 41 .4 s 41 .3 s 41 .2 s 41 .1 s Fig. 1 | The magnetic field orientation of the gravitationally lensed galaxy 9io9 at z = 2.553. a–d, ALMA 242 GHz polarimetric observations of the Stokes I, Q and U parameters, and the polarized intensity (PI). The synthetic beam of the observations (1.2″ × 0.9″, θ = 68°) is shown as the red ellipse, lower left. The B field orientation is indicated by white lines shown at the Nyquist sampling, with line lengths proportional to the polarization fraction. e–h, Synthetic polarimetric observations using a constant B field configuration in the source plane. Contours indicate signal to noise: for Stokes I, the contours increase as σI × 2 3,4,5,…. For Stokes Q and U and for PI, the contours start at 3σ and increase in steps of 1σ. Dec., declination; RA, right ascension. a b cSource plane Lensed source Synthetic observation –0.4 –0.2 0 0.2 0.4 Offset (″) –0.4 –0.2 0 0.2 0.4 O ff se t (″ ) Psp = 1.0%FB,sp = 5º –4 –2 0 2 4 –4 –2 0 2 4 –4 –2 0 2 4 –4 –2 0 2 4 O ff se t (″ ) Offset (″) Offset (″) O ff se t (″ ) 〈PI〉 = 1.0 ± 0.0% 〈FB,I〉 = 5.0 ± 0.0º 〈Pip〉 = 1.1 ± 0.1%〈FB,ip〉 = 4.3 ± 6.3º Fig. 2 | Source plane configuration of the magnetic field and lensing model. a, Source plane intensity and field orientation. b, Lensed source plane image. c, Synthetic observations with the synthetic beam size (1.2″ × 0.9″, θ = 68°) indicated by the red ellipse. The B field orientation is indicated by white lines with lengths proportional to the polarization fraction. The median and root mean squared values of the polarization fraction and B field orientation are indicated in the images. The caustics in the source plane and image plane are shown as green and yellow lines, respectively. Nature | www.nature.com | 3 and the star-formation rate density exceeds 100 M⊙ yr −1 kpc−2 (ref. 5). The high dense gas fraction of the molecular reservoir—as traced by the ratio of CO(4–3)/C i(1–0) emission—is also consistent with the injection of supersonic turbulence, which plays a key role in shaping the lognormal probability distribution function of the molecular gas density20. There is also tentative evidence of stellar feedback in action through the broad lines of dense gas tracers5. Finally, one expects a high cosmic-ray flux density in the ISM of 9io9, commensurate with the high star-formation rate density, and this too could serve to amplify magnetic fields. Therefore, 9io9 probably has the conditions required to rapidly amplify any weak seed fields by means of the small-scale dynamo effect, with amplification occurring on scales up to and includ- ing the full star-forming disk. Assuming equipartition between the turbulent kinetic and magnetic energies, we estimate an upper limit of the equipartition turbulent B field strength of 514 μG (Methods). This is comparable to the estimated turbulent B field strength of 305 ± 15 μG within the central kiloparsec of the starburst region of M82 also using far-infrared polarimetric observations21. This indicates that the star- burst activity of 9io9 could be be driving the amplification of B fields across the disk. Feedback-induced turbulence is a route to accelerating the growth of the seed fields, but to produce the ordered field on the kpc-scales observed requires a mean-field dynamo14,22. This mean-field dynamo can be achieved through the rapid differential rotation of the gas disk, and this provides a mechanism for the ordering of an amplified B field driven by star formation and stellar feedback processes. 9io9 is turbulent, intensely star-forming and rapidly rotating (maximum velocity vmax ≅ 300 km s−1). This suggests that rather than an episode of violent feedback priming a large-scale but turbulent field that later evolves into an ordered field during a period of relative quies- cence19, the small-scale and mean-field dynamo mechanisms oper- ate in tandem. We estimate that the mean-field dynamo in 9io9 has not yet had time to maintain or amplify the B field (Methods). This indicates that the intense starburst is most important in amplify- ing the galactic field at z = 2.6. We postulate that this ‘dual dynamo’ might be the common mode by which galactic-scale ordered mag- netic fields are established in young gas-rich, turbulent galaxies in the early Universe. Coherent magnetic fields consistent with the mean-field dynamo have been observed at z = 0.4 by means of Faraday rotation of a back- ground polarized radio source3 (note that such observations are not possible for 9io9). Magnetic fields are already known to be present in the environment around normal galaxies at z ≅ 1 as revealed by the association of Mg ii absorption systems along quasar sightlines that show Faraday rotation23, and indirectly through the existence of radio synchrotron emission from star-forming galaxies. However, mapping the B fields in individual galaxies at high redshift has, so far, proved challenging. Our observations show that the polarized emission from magnetically aligned dust grains is a powerful tool to trace the B fields of the cold and dense ISM in high-redshift galaxies. Galaxy 9io9 is a particularly luminous example of a population of dusty star-forming galaxies in the early Universe that contribute to the cosmic infrared background (CIB). If the 1% level of polarization detected in 9io9 is representative of the general population of dusty star-forming galaxies24 then routine detection and mapping of mag- netic fields in galaxies at high redshift is feasible (that is, in integration times of less than 24 h) even in unlensed systems with ALMA. This offers a new window to characterize the physical conditions of the ISM in galaxies when galaxy growth was at its maximum, and will enable a better understanding the role of magnetic fields in shaping the early stages of galaxy evolution. The strength of the galactic magnetic field in local spiral galaxies is of order 10 μG (ref. 1), and up to an order of magnitude higher in starbursts8. Without resolving the polarization field in 9io9 below 100 pc scales it is not possible to reliably estimate the B field strength using dust polarization observations. Nevertheless, given the injection of kinetic turbulence driven by stellar feedback we estimate the strength of the B field in 9io9 to be probably greater than that of local spiral galaxies, but similar to that of the central regions of nearby starburst galaxies (Methods). Finally, these observations imply the CIB itself may be weakly polar- ized24,25. Although misalignments of galaxies along the line of sight will serve to reduce the net polarization of the CIB, if the orientation of disks that host galactic-scale ordered B fields is correlated on large scales due to tidal alignments26, then a polarization signal could remain and therefore fluctuations in the polarization intensity of the CIB could be used as a new probe of the physics of structure formation25. This has consequences for cosmological experiments that seek to derive infor- mation on primordial conditions from observations of the polarization of the cosmic microwave background, especially if a curl component is present in the CIB polarization field25,27. A polarized component of the CIB at millimetre wavelengths, of extragalactic origin, dominated by emission at z ≅ 2 and with a power spectrum that is shaped by large-scale structure at this epoch, will be a subtle but important foreground for future precision cosmic microwave background experiments to con- tend with. Online content Any methods, additional references, Nature Portfolio reporting summa- ries, source data, extended data, supplementary information, acknowl- edgements, peer review information; details of author contributions and competing interests; and statements of data and code availability are available at https://doi.org/10.1038/s41586-023-06346-4. 1. Beck, R., Chamandy, L., Elson, E. & Blackman, E. G. Synthesizing observations and theory to understand galactic magnetic fields: progress and challenges. Galaxies 8, 4 (2019). 2. Lopez-Rodriguez, E. et al. Extragalactic magnetism with SOFIA (SALSA Legacy Program): the magnetic fields in the multiphase interstellar medium of the antennae galaxies. Astrophys. J. Lett. 942, L13 (2022). 3. Mao, S. A. et al. Detection of microgauss coherent magnetic fields in a galaxy five billion years ago. Nat. Astron. 1, 621–626 (2017). 4. Geach, J. E. et al. The Red Radio Ring: a gravitationally lensed hyperluminous infrared radio galaxy at z = 2.553 discovered through the citizen science project Space Warps. Monthly Notices R. Astronom. Soc. 452, 502–510 (2015). 5. Geach, J. E., Ivison, R. J., Dye, S. & Oteo, I. A magnified view of circumnuclear star formation and feedback around an active galactic nucleus at z = 2.6. Astrophys. J. Lett. 866, L12 (2018). 6. Hoang, T. & Lazarian, A. Radiative torque alignment: essential physical processes. Mon. Not. R. Astronom. Soc. 388, 117–143 (2008). 7. Andersson, B.-G., Lazarian, A. & Vaillancourt, J. E. Interstellar dust grain alignment. Ann. Rev. Astron. Astrophys. 53, 501–539 (2015). 8. Lopez-Rodriguez, E. et al. Extragalactic magnetism with SOFIA (SALSA Legacy Program). IV. Program overview and first results on the polarization fraction. Astrophys. J. 936, 92 (2022). 9. Saintonge, A. & Catinella, B. The cold interstellar medium of galaxies in the local universe. Annu. Rev. Astron. Astrophys. 60, 319–361 (2022). 10. Beck, R., Poezd, A. D., Shukurov, A. & Sokoloff, D. D. Dynamos in evolving galaxies. Astron. Astrophys. 289, 94–100 (1994). 11. Brandenburg, A. & Subramanian, K. Astrophysical magnetic fields and nonlinear dynamo theory. Phys. Rep. 417, 1–209 (2005). 12. Schober, J., Schleicher, D. R. G. & Klessen, R. S. Magnetic field amplification in young galaxies. Astron. Astrophys. 560, A87 (2013). 13. Biermann, L. Über den Ursprung der Magnetfelder auf Sternen und im interstellaren Raum (miteinem Anhang von A. Schlüter). Zeitschrift Naturforschung Teil A 5, 65 (1950). 14. Beck, R., Brandenburg, A., Moss, D., Shukurov, A. & Sokoloff, D. Galactic magnetism: recent developments and perspectives. ARA&A 34, 155–206 (1996). 15. Kulsrud, R. M. & Zweibel, E. G. On the origin of cosmic magnetic fields. Rep. Progr. Phys. 71, 046901 (2008). 16. Subramanian, K. From primordial seed magnetic fields to the galactic dynamo. Galaxies 7, 47 (2019). 17. Lee, H. M. & Draine, B. T. Infrared extinction and polarization due to partially aligned spheroidal grains: models for the dust toward the BN object. Astrophys. J. 290, 211–228 (1985). 18. Rieder, M. & Teyssier, R. A small-scale dynamo in feedback-dominated galaxies as the origin of cosmic magnetic fields - I. The kinematic phase. Mon. Not. R. Astron. Soc. 457, 1722–1738 (2016). 19. Rieder, M. & Teyssier, R. A small-scale dynamo in feedback-dominated galaxies—III. Cosmological simulations. Mon. Not. R. Astron. Soc. 472, 4368–4373 (2017). 20. Padoan, P., & Nordlund, A. The stellar initial mass function from turbulent fragmentation. Astrophys. J. 576, 870–879 (2002). 4 | Nature | www.nature.com Article 21. Lopez-Rodriguez, E., Guerra, J. A., Asgari-Targhi, M. & Schmelz, J. T. The strength and structure of the magnetic field in the galactic outflow of Messier 82. Astrophys. J. 914, 24 (2021). 22. Ruzmaikin, A., Sokolov, D. & Shukurov, A. Magnetism of spiral galaxies. Nature 336, 341–347 (1988). 23. Bernet, M. L., Miniati, F., Lilly, S. J., Kronberg, P. P. & Dessauges-Zavadsky, M. Strong magnetic fields in normal galaxies at high redshift. Nature 454, 302–304 (2008). 24. Bonavera, L., González-Nuevo, J., De Marco, B., Argüeso, F. & Toffolatti, L. Statistics of the fractional polarization of extragalactic dusty sources in Planck HFI maps. Mon. Not. R. Astron. Soc. 472, 628–635 (2017). 25. Feng, C. & Holder, G. Polarization of the cosmic infrared background fluctuations. Astrophys. J. 897, 140 (2020). 26. Hirata, C. M. & Seljak, U. Intrinsic alignment-lensing interference as a contaminant of cosmic shear. Phys. Rev. D 70, 063526 (2004). 27. Lagache, G., Béthermin, M., Montier, L., Serra, P. & Tucci, M. Impact of polarised extragalactic sources on the measurement of CMB B-mode anisotropies. Astron. Astrophys. 642, A232 (2020). Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. © The Author(s) 2023 Methods Observations and data reduction Galaxy 9io9 at celestial coordinates α = 02 h 09 min 41.3 s, δ = 00° 15′ 58.5″ ( J2000) was observed by the ALMA 12 m array in project 2021.1.01461.S. The target was observed in two sessions, in which each session consisted of observations of sufficient duration to measure the rotation of the parallactic angle of the telescope by 60° or more. The first session, comprising two Execution Blocks, was executed on 12 April 2022, and the second, comprising one Execution Block, was executed on 14 April 2022. The ALMA grid source J2253 + 1608 was observed as a bandpass and flux calibrator, J0006–0623 was used as the polarization calibrator, J0208–0047 was used as the phase calibrator for the first session and J0217 + 0144 was used as the phase calibrator for the second session. The data were processed in the Common Analysis Software Applica- tions (CASA) v.6.2.1 using both the standard ALMA calibration pipeline and an extra polarization calibration script. The standard calibration pipeline first applied a series of steps that included corrections to the amplitudes (based on system temperature measurements) and phases (based on water vapour radiometer measurements). This was followed by calibration steps in which corrections for the amplitudes and phases as a function of frequency were derived using the bandpass calibrator, scaling factors for the amplitudes were derived using the flux calibra- tor and corrections for the phases and amplitudes versus time were derived using the phase calibrator. The polarization calibration script derived and applied further corrections for instrument-related polari- zation effects related to imperfections in the feeds and the orientation of the feeds as a function of time. The uncertainty in the final linear polarization calibration fraction was 0.1%, whereas the uncertainty in the polarization angle was 1°. We also used CASA for imaging of the calibrated measurement set. The procedure tclean was used to produce the IQUV Stokes images, with the ‘clarkstokes’ deconvolver, natural weighting of the visibilities and a common restoring beam (ensuring all Stokes images share the same synthetic beam). Dirty images were cleaned to a stopping threshold of 10 μJy and the final beam shape had a full-width at half- maximum of 1.2″ × 0.9″ at position angle 68°. All images were primary beam-corrected to account for fall-off in sensitivity away from the phase centre, but note that the target was compact and this represented a negligible correction to the measured flux densities. ALMA provided a systematic uncertainty in linear polarization of 0.03% with minimum detectable polarization of 0.1%, and therefore all quoted uncertainties in this work have the minimum detectable polarization of 0.1% added in quadrature. Measuring polarization To account for the vector quantity of the polarization measurements, we estimate the integrated polarization fraction as P Q U b I =   +   −   (1) 2 2 and the B field orientation as χ U Q = 1 2 arctan     + pi 2 (2)B       where ⟨I⟩, ⟨Q⟩ and ⟨U⟩ are the mean of the Stokes I, Q and U for measure- ments with a polarized intensity signal to noise ratio P/σP ≥ 3, and bias b is the standard error of the uncertainties of the Stokes Q and U images, σQ and σU. The mean polarization fraction, ⟨P⟩, and B field orientation, ⟨χB⟩, are estimated as the mean distribution of the individual measure- ments of P and χB in independent beams (that is, using the individual Stokes Q and U and then computing P and χB). The uncertainties in these values are estimated as the standard deviation of each distribution, respectively. Extended Data Fig. 1 shows the distribution of individual measurements at the Nyquist sampling for pixels with PI/σPI ≥ 3. Lens modelling Our goal is to obtain the configuration of the B field morphology in the source plane using previous knowledge of the lens model of 9io9. The lens model is computed using a singular isothermal ellipsoid and a shear component. The lens is described by the critical radius (that is, Einstein radius), θE, lens offset position from the central coordinates, xc, yc, ellipticity, ϵ, shear, γ, and shear orientation, θγ. Extended Data Table 1 summarizes the best-fit lens model parameters from ref. 5. The morphology of the continuum thermal emission in the source plane is assumed to match that of the CO(4–3) molecular gas emission5. The CO(4–3) spectral cube can be fit with a kinematic model based on the ring or disk morphologies of the molecular gas reservoirs of local ultraluminous galaxies28, with the source plane spectral cube well modelled with a disk of maximum radius R = 322max −20 +11  mas (2, 647−160 +88  pc) and maximum rotation speed of V = 300max  km s −1). The disk has a tilt angle of 5 ± 4° (east of north) and is inclined along the line of sight by 50−8 +3°. We produce synthetic polarimetric observations of our ALMA pola- rimetric observations as follows. First, we compute the singular iso- thermal ellipsoid and shear lens model from the parameters given in Extended Data Table  1 using LENSTRONOMY29. Using a grid of 801 × 801 pixels2 with a scale of 0.01″ pixel−1, which corresponds to a spatial scale of 82 pc in the source plane, we model the source plane using an asymmetric 2D Gaussian profile with a full-width at half-maximum equal to the Rmax and Rmin at an angle of θ (Extended Data Table 1). The source plane and image plane at the native resolution of 0.01″ is shown in Fig. 2, also showing the caustics of the lens model. To mimic the observed data, we use the simobserve task of CASA, which simulates the observation of a given sky model (that is, the modelled image plane) with the ALMA 12 m array. We match the parameters of the real observations as closely as possible: we ‘observe’ in two sessions using antenna configuration ‘2’ for Cycle 8. The integrations start at hour angles of −2.339 and −2.422 h for the two sessions and last 1.86 and 0.86 h on 12 and 14 April 2022, respectively, ensuring an identical sampling of the uv plane. A model is applied to simulate noise in the simulated observations, dominated by a thermal component, using the atmospheric transmission at microwaves model to simulate the atmospheric profile at the ALMA site, with the precipitable water vapour column as a scaling parameter. We assume a precipitable water vapour of 0.6 mm, equivalent to the average column across the obser- vations. We clean the simulated visibilities in the same manner as the real data, but apply a restoring beam equivalent to the common beam derived for the real observations, to ensure that the final synthesized beams of the simulated images match the real data. The B field orientation in the source plane (sp) is assumed to be con- stant, χB,sp, with a constant polarization fraction of unity, Psp. The model Stokes Qsp and Usp in the source plane are computed as Q P χ I= cos(2 ) × (3)Bsp sp ,sp sp U P χ I= sin(2 ) × (4)Bsp sp ,sp sp where Isp is the Stokes I in the source plane. We multiply by Isp to convert the Stokes Qsp and Usp in surface brightness, which allows us to compute the lens model. Figure 2 shows the B field orientation in the source plane with vector line lengths proportional to the polarization fraction, that is, Psp = 1%. As the Stokes Qsp and Usp are density profiles, these images can also be lensed using the same procedure as the Stokes Isp. The polarized intensity, PIip, polarization fraction, Pip, and B field orientation, χB,ip, of the final synthetic polarimetric observation in the image plane are computed as Article Q UPI = + (5)ip ip 2 ip 2 P I = PI (6)ip ip ip      χ U Q = 1 2 arctan + pi 2 (7)B,ip ip ip The final synthetic observations Stokes Iip, Qip and Uip are shown in Fig. 1. This figure also shows the synthetic polarized intensity, PIip, and the B field orientation, χB,ip, with the length of the lines proportional to the synthetic polarization fraction Pip. Note that the B field orientation of the synthetic observations does not have the added rotation of 90° as shown in equation (2). This is because we are modelling the B field rather than the E vector measured by the ALMA polarimetric observa- tions. However, to compare with our observations, Fig. 1 shows the synthetic Stokes Q and U in the same reference (that is, E vector) as the ALMA polarimetric observations. Constraining the B field orientation It has been shown that for a non-rotating lens the polarization vector of the electromagnetic wave does not rotate. This result has been dem- onstrated using pure geometrical definitions30,31, and from first princi- ples using a Newtonian potential and solving the spacetime metrics32. In general, the photons travel along null geodesics, which indicates that the vector properties are time invariant across the geodesic. This property allows us to study the intrinsic B field geometry of the source amplified by the gravitational lens. We explore potential changes of the B field orientation in our lens model. With a polarization fraction of 1% we set a range of constant χB,sp orientations in the source plane. Extended Data Fig. 2 shows the source plane and lensed model for several values of χB,sp = 5, 10, 45, 90° at a constant Psp = 1% using the same lens model as described above, and simulating the observations in the same way to produce convolved and noisy images at the same scale as the data. We conclude that the lens model and convolution do not produce rotation or change in the polarization fraction from the source plane to the image plane. For χB,sp = [5, 10]°, the final Bip has an orientation similar to the beam position angle in regions with a low signal to noise ratio. However, the median B field orientation in the final synthetic observation is consistent with the modelled B field orientation in the source plane, with the presence of noise and asymmetric synthetic beam contrib- uting to the angular dispersion in the image plane. We find that the background noise and uv plane sampling (in general, imaging of the visibilities) introduce an uncertainty of roughly 10° in ⟨χB,ip⟩, in which the uncertainty is dominated by the signal to noise ratio of the ALMA polarimetric observations. The B field orientation can also be constrained using the Stokes Q and U images. Our observations show that Stokes Q is negative and Stokes U is consistent with zero (Fig. 1). This result shows that the E vector is mainly in the north–south direction. This configuration gives us the opportunity to tightly constrain the possible range of values in the B field orientation. Extended Data Fig. 3 shows the Stokes I, Q and U images, and the polarized intensity of synthetic observations for the same set of artificial B field orientations in the source plane. Note the change from negative to positive Q from χB,sp of 5 to 90°. This behaviour is expected, but we show it for completeness. We can constrain the orientation by rotating the B field until the corresponding Stokes U image shows a 3σ detection in the synthetic polarimetric observations. This condition is met when χB,sp = ± 10°. For orientations deviating out- side ±10° the synthetic observations are not consistent with our obser- vations. We conclude that the most likely B field orientation consistent with our observations is χ = 5B,sp −15 +5°, where we have set the B field ori- entation to be parallel to the CO(4–3) emission of the disk in the source plane5. Constraining the polarization fraction We assume a constant polarization fraction of 1% in the source plane. The lens model shows that the median polarization fraction is con- sistent with this level without dispersion (Extended Data Fig. 2). We conclude that the lens model does not change the polarization frac- tion from the source plane to the image plane. Using the simulated observations, we estimate the median and r.m.s. polarization fraction, finding that the combination of uv-plane sampling and background noise adds an uncertainty of roughly 0.1–0.2% in the polarization fraction. Energy equipartition We estimate the B field strength assuming equipartition between the turbulent kinetic energy and the turbulent magnetic energy. Let the turbulent kinetic energy, UK, and turbulent magnetic energy, UB, be U ρσ= 1 2 (8)K v 2 U B pi = 8 (9)B 2 where ρ is the volume density, σv is the velocity dispersion and B is the field strength. Then, assuming equipartition, UK = UB, the field strength is B ρ σ= 4pi . (10)eq v To estimate the baryonic volume density, ρ = M/V, we use the disk volume, V = πr2h, and mass of the molecular gas of the galaxy. The radius is r = 2, 647−160 +88  pc and the height is h ≤ 600 pc (ref. 5); note that the disk height is unresolved in the observations of the CO molecular gas, so we take the size of the beam of the observations as a conserva- tive upper limit, also noting that the scale heights of bursty disk galax- ies at high redshift are observed to be larger than local disks, with heights of several hundred parsecs typical. The total molecular mass is estimated from the CO luminosity5 as M = (7.5 ± 0.1) × 10H 10 2  M⊙ and the velocity dispersion of the gas is σv = 73 ± 4 km s −1, from kinematic fitting of the disk model to the CO data cube5. From this, we estimate an equipartition B field strength of Beq ≤ 514 μG, however note that the observed velocity dispersion may be affected by large-scale flows from galactic winds and shearing effects by the rotating galactic disk, such that the turbulent component velocity dispersion might be lower. These effects produce an overestimated measurement of the B field strength33. Spiral galaxies have an average ordered B field strength of around 5 ± 2 μG with a total B field strength of 17 ± 14 μG assuming equiparti- tion between the total B field and total cosmic-ray electron density1. A revised equipartition formula to account for energy losses in nearby (within 160 Mpc) starburst galaxies estimated34 equipartition B field strengths in the range 70–770 μG. Recent far-infrared polarimetric observations21 of the starburst galaxy M82 showed that the turbulent B field strength is 305 ± 15 μG within the central kiloparsec region. These authors modified the Davis–Chandrasekhar–Fermi method to account for the large-scale flow of the galactic outflow. The correction from the equipartition B field strength, Beq = 540 ± 170 μG, was estimated to be roughly 25%. Ordering timescale We estimate the timescale to order a large-scale magnetic field in a disk. The ordering timescale can be estimated as35 ∣ ∣∼t h l βh D= (11) 2 c d −1/4 where h is the half-thickness of the galactic disk and lc is the large-scale coherence length. β is the turbulent diffusivity defined as β = lσv/3, where l is the coherence length of the small-scale turbulence and σv is the turbulent velocity. Dd is the dynamo number defined as      ∣ ∣D hΩ σ = 9 (12)d v 2 where Ω is the angular velocity. For 9io9 we assume a Gaussian scale height5 h = 300 pc and σv = 73 ± 4 km s −1 of the molecular gas. To estimate the angular velocity we use the deprojected circular velocity, v = 360max −11 +49  km s−1, at the maximum radius of the disk, r = 2, 647max −160 +88  pc, yielding Ω v r= / = 145 ± 8max max km s −1 kpc−1. The typical coherence length of the small-scale dynamo is driven by stellar activity in the molecular gas of galaxies on scales of l = 1–10 pc (ref. 36). We assume a large-scale coherence length in the range lc = 0.5 − 2.0 kpc, corresponding to 50–2,000 times larger than the turbulence coherence length, l. We estimate a dynamo number Dd = 3.2 ± 0.4 and turbulent diffusivity β = (4.2 ± 2.0) × 10 25 cm2 s−1. Finally, the ordering timescale is estimated to be in the range of t ≅ 0.4–18 Gyr, noting that only 2.5 Gyr have elapsed by z = 2.6. Analyti- cal solutions of the evolution of B fields in spiral galaxies predicts that a large-scale B field can be formed within 2–5 Gyr (ref. 37). This study suggests that the mean-field dynamo may be already active at z < 3. Note that 9io9 has a very high star-formation rate, which may inject higher turbulent energy into the system than those studied in ref. 37 affecting the ordering of the B field in the galaxy’s disk. The regular B field is amplified if the dynamo number is larger than the critical value of Dd,cr ≅ 7 estimated from numerical simulations of galactic dynamo models22. For comparison, the Milky Way has a dynamo number of Dd ≅ 9 > Dd,cr, which shows that large-scale dynamo mechanism is important in our galaxy. For 9io9 we estimate Dd = 3.2 ± 0.4, below the critical value. In addition, the time for 9io9 to complete a rotation is t ≅ 45 Myr, which indicates that roughly 9–400 galactic rotations at a galactocentric radius of 2.6 kpc are required to order the large- scale B field. Data availability The ALMA data used in this work will be available through the ALMA Science Archive (https://almascience.nrao.edu/aq) with reference to the project code 2021.1.01461.S. 28. Downes, D. & Solomon, P. M. Rotating nuclear rings and extreme starbursts in ultraluminous galaxies. Astrophys. J. 507, 615–654 (1998). 29. Birrer, S. & Amara, A. lenstronomy: multi-purpose gravitational lens modelling software package. Phys. Dark Universe 22, 189–201 (2018). 30. Dyer, C. C. & Shaver, E. G. On the rotation of polarization by a gravitational lens. Astrophys. J. Lett. 390, L5 (1992). 31. Burns, C. R. Gravitational Lensing of Polarized Sources. PhD thesis, Univ. Toronto (2002). 32. Faraoni, V. On the rotation of polarization by a gravitational lens. Astron. Astrophys. 272, 385 (1993). 33. Guerra, J. A., Lopez-Rodriguez, E., Chuss, D. T., Butterfield, N. O. & Schmelz, J. T. The strength of the sheared magnetic field in the Galactic’s circumnuclear disk. Astron. J. 166, 37 (2023). 34. Lacki, B. C. & Beck, R. The equipartition magnetic field formula in starburst galaxies: accounting for pionic secondaries and strong energy losses. Mon. Not. R. Astron. Soc. 430, 3171–3186 (2013). 35. Arshakian, T. G., Beck, R., Krause, M. & Sokoloff, D. Evolution of magnetic fields in galaxies and future observational tests with the Square Kilometre Array. Astron. Astrophys. 494, 21–32 (2009). 36. Haverkorn, M., Brown, J. C., Gaensler, B. M. & McClure-Griffiths, N. M. The outer scale of turbulence in the magnetoionized galactic interstellar medium. Astrophys. J. 680, 362–370 (2008). 37. Rodrigues, L. F. S., Chamandy, L., Shukurov, A., Baugh, C. M. & Taylor, A. R. Evolution of galactic magnetic fields. Mon. Not. R. Astron. Soc. 483, 2424–2440 (2019). Acknowledgements This paper makes use of the following ALMA data: ADS/JAO. ALMA#021.1.01461.S. ALMA is a partnership of the European Southern Observatory (ESO) (representing its member states), National Science Foundation (USA) and National Institute of Natural Sciences (Japan), together with the National Research Council (Canada), Ministry of Science and Technology (MOST) and Academia Sinica Institute of Astronomy and Astrophysics (ASIAA) (Taiwan), and Korea Astronomy and Space Science Institute (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, Associated Universities Inc. (AUI)/National Radio Astronomy Observatory (NRAO) and National Astronomical Observatory of Japan (NAOJ). J.E.G. and M.J.D. acknowledge support from the Royal Society. E.L.-R. is funded by the Universities Space Research Association, Inc. under the 07_0032 NASA/DLR Stratospheric Observatory for Infrared Astronomy programme. G.J.B. acknowledges support from Science and Technology Facilities Council grant no. ST/T001488/1. Author contributions J.E.G. led the proposal to obtain the data and designed the observations. J.E.G. and G.J.B. reduced the data. J.E.G. and E.L.-R. performed the analysis. All authors contributed to the manuscript. Competing interests The authors declare no competing interests. Additional information Correspondence and requests for materials should be addressed to J. E. Geach. Peer review information Nature thanks the anonymous reviewers for their contribution to the peer review of this work. Reprints and permissions information is available at http://www.nature.com/reprints. Article Extended Data Fig. 1 | Histograms of the measured polarization fraction and magnetic field orientation. The polarization fraction (a) is displayed in bins of 0.5% and the B field orientation (b) in bins of 10 degrees, where individual measurements are taken in independent beams. The median and 1σ scatter of the distributions is given for each measurement. Extended Data Fig. 2 | Characterisation of the magnetic field orientation and polarization fraction. The B field orientation (white lines) at 5 degrees (a), 10 degrees (b), 45 degrees (c), and 90 degrees (d) are shown over the source plane (a), lens model (e), and synthetic observations (f). We indicate the caustics in the source plane (green line), and image plane (yellow line). The median and r.m.s. of the polarization fraction and B field orientation are shown at the top right of each panel. The red ellipse indicates the synthetic beam of the simulated observations. Article Extended Data Fig. 3 | Variation of the Stokes Q and U parameters with magnetic field orientation. The B field orientation (white lines) at 5 degrees (a), 10 degrees (b), 45 degrees (c), 90 degrees (d) are shown over the Stokes I (a), Stokes Q (e), Stokes U (f), and polarized intensity (g) of the synthetic polarimetric observations. The red ellipse shows the  synthetic beam of the simulated observations (a). For Stokes I, the contours increase as σI × 2 3,4,5,…. For Stokes Q and U, and PI, the contours start at 3σ and increase in steps of 1σ. Note that at σB,sp = 10 degrees a 3σ signal becomes detectable in Stokes U. Extended Data Table 1 | Lens model parameters Angles are measured counterclockwise, i.e., East of North. Coordinates in J2000.