Young Radio Sources Expanding in Gas-Rich ISM: Using Cold Molecular Gas to Trace Their Impact

We present the results from the study of the resolved distribution of cold molecular gas around eight young (<10^6 yr), peaked-spectrum radio galaxies. This has allowed us to trace the interplay between the radio jets and the surrounding medium. For three of these sources we present new CO(1-0) observations, obtained with NOEMA. In two targets, we also detected CN lines, both in emission and absorption. Combining the new observations with already published data, we discuss the main results obtained. Although we found that a large fraction of the cold molecular gas is distributed in disc-like rotating structures, in most of the sources high turbulence and deviations from purely quiescent gas (including outflows) were observed in the region co-spatial with the radio continuum emission. This suggests the presence of an interaction between radio plasma and cold molecular gas. We found that newly born and young radio jets, even those with low power i.e., P_jet<10^45 erg/s), can drive massive outflows of cold, molecular gas. The outflows are, however, limited to the sub-kpc regions and likely short lived. On larger scales (a few kpc), we observed cases where the molecular gas appears to avoid the radio lobes and, instead, wraps around them. The results suggest the presence of an evolutionary sequence, consistent with simulations, where the type of impact of the radio plasma changes as the jet expands, going from a direct jet-cloud interaction on sub-kpc scales to a gentler pushing aside of the gas, increasing its turbulence and likely limiting its cooling. This effect can be mediated by the cocoon of shocked gas inflated by the jet-cloud interactions. Building larger samples of young and evolved radio sources for observation at a similar depth and spatial resolution to test this scenario is now needed and may be possible thanks to more data becoming available in the growing public archives.


Introduction
Tracing the properties of the gas in the centre of galaxies has provided important insights into the interplay between the energy emitted by an active galactic nucleus (AGN) and the surrounding medium. Two highly intertwined processes, fuelling and feedback, play a main role (e.g., [1,2]) and are thought to be crucial in the evolution of the host galaxy. Among the various types of AGNs, radio (jetted) AGNs are known to have jets and lobes ranging in size from pc to hundreds of kpc. Spanning such large scales implies that this type of AGN can arXiv:2302.14095v1 [astro-ph.GA] 27 Feb 2023 potentially have an impact on the gas in the interstellar medium (ISM), circum-galactic medium (CGM), and intergalactic medium (IGM). They are, therefore, particularly interesting for the study of AGN feedback if the radio plasma efficiently couples with this gas. Although the effect of jets and lobes on larger, intergalactic scales is well established, thanks to studies of X-ray cavities ( [3,4]), it is now becoming clear that radio jets can also impact the ISM inside the host galaxy. As suggested by numerical simulations (e.g., [5]), this can occur in different ways, from an initial phase with a strong interaction (possibly driving fast and massive outflows) to a more rapid, momentum-dominated expansion phase when the jet breaks free from the denser central regions and reaches galactic sizes and beyond.
The actual impact will, however, depend on several parameters: the jet power, orientation of the jet with respect to the gaseous structure, and evolutionary stage of the jet. This complexity is well described by the results from recent numerical simulations (e.g., [5][6][7][8][9][10]). In particular, the simulations developed by [5,10] and more recently expanded by [7] illustrate the process of a jet entering a clumpy medium and initially directly interacting with the surrounding gas clouds. They also show how the impact of the jet is amplified by the creation of a cocoon of a shocked medium resulting from the jet-cloud interactions and expanding orthogonal to the jet (see Figure 1 for an example). This latter phase starts in the galactic ISM and is also visible as an X-ray emission around the radio source (e.g., [11][12][13][14]), bearing similarities with the X-ray cavities observed in the IGM.
Some numerical simulations are also starting to provide direct comparisons to observations [15,16]. An example illustrating how an interaction proceeds is shown in Figure 1, where simulations recently proposed to describe the results from the cold molecular gas for the young radio galaxy B2 0258+35 are shown [16] (see also Section 4 for details).   [7], showing the evolution of the jet-disc system in the case of a jet slightly offset compared to the disc. The figure is taken from [16] and describes the case of the radio galaxy B2 0258+35. The jet plasma is shown in blue and the dense ISM clouds in orange-red. The ablated gas and shocked ambient medium are in yellow. The strongest interactions occur within the disc where the main jet hits the clouds head-on. These regions show enhanced velocity dispersion and bulk velocities of a few hundred km s −1 and are locations of sharp jet deflection and splitting. The outer disc is dispersed but remains largely intact (see [7,16]).
Thus, both the simulations and several observational results suggest that the initial phase (<10 6 yr; [17][18][19]) of the evolution of the radio jet, with the jet expanding through the inner (few) kpc of the host galaxy, can be particularly important. Tracing the gas around the young radio jets is, therefore, key to quantifying the effects of this kind of interaction, which were described above. This has been carried out extensively for the warm ionised gas ( [20][21][22]). However, the cold phase of the gas -and, in particular, the cold molecular gas that is used in this paper -is highly relevant for its direct connection to star formation. Furthermore, it can be observed with the high spatial resolution necessary for tracing this interaction.
Here, we provide an overview of the results obtained so far by an ongoing project aimed at tracing the properties of the cold gas in young radio galaxies embedded in a rich gaseous medium. We focus on the properties of the cold molecular gas component traced with sensitive, high spatial resolution observations carried out with the Northern Extended Millimeter Array (NOEMA) and the Atacama Large Millimeter Array (ALMA). The high spatial resolution obtained for most of the targets has been used to resolve the distribution and kinematics of the gas across the radio emission and to better describe the effects (if any) of the jet-ISM interaction.
This paper is structured as follows. In Section 3, we describe the selection of the target group of objects of our observations and our study; in Section 3, we present the as-yet unpublished CO(1-0) results for the three objects; in Section 4, we summarise the published results for the rest of the sample; in Section 5, we combine our findings and highlight the main points, and in Section 6, we finish with some conclusions about the implications of the results obtained so far and remarks about future work.
Throughout this paper we adopt a cosmology that assumes a flat Universe and the following parameters:

Cold Molecular Gas in Young Radio Sources: Overview of Our Observations
The eight sources observed in CO are listed in Table 1, together with a selection of their relevant properties. The list includes known young radio galaxies (i.e., small and peaked spectrum sources) residing in a gas-rich medium. In particular, the selected objects have all been detected in H I absorption, ensuring the presence of cold gas [16,[24][25][26][27][28].
The targets are all small radio sources of linear size 4 kpc, with a peaked radio spectrum and a turnover frequency ranging from about 100 MHz to a few GHz. These are properties typical of radio galaxies recently born (i.e., 10 6 yr old; see [18] for a review). In all cases, VLBI observations of their radio continuum structures are available in the literature.
The original goal of the project was to build a sample covering a broad range of source ages and jet powers (using the radio power as a proxy to a first approximation) to be able to investigate the trends proposed by the numerical simulations. Although the final list of observed objects is still limited, the observed objects cover a relatively large range in these parameters. In more detail, the ages cover a range from 91 yr (obtained from hotspot advance measurements [29]) in the case of PKS 1718-63, which is considered one of the youngest known radio sources, to 3.5 Myr for 3C 305 [30], the oldest in the sample. Additionally, the radio power (proxy for jet power) of the sources in our sample covers a large range, as seen in Table 1, with PKS 1718-63, IC 5063, and B2 0258+35 at the lower end of the distribution and PKS 0023-26 and PKS 1549-79 the most powerful sources in the sample. These powerful sources also host powerful radiative AGNs [20,21]. The presence of both radiative and mechanical energy from these AGNs makes them good candidates for tracing powerful feedback effects. For more details, we also refer to the papers listed in Table 1.
Of the eight sources, in five of them, the cold molecular gas has been extensively discussed (see references in Table 1 for details and Section 4 for a short overview of recent results). For the three remaining sources, in Section 3, we present the results from new NOEMA CO(1-0) observations. For one object (4C 31.04), a paper is currently being prepared that includes, together with CO, new H I VLBI observations (Murthy et al. in prep). As a result, here, we present only the highlights of the CO results relevant to the discussion in this paper. The columns are: (1) source name; (2) telescopes used; (3) beam of the CO observations (in parenthesis, the scale in kpc); (4) the radio luminosity at 1.4 GHz; (5) the size of the radio sources (in the case of B2 0258+35, only the central, young source is given, see [23] for details); (6) the structure of the distribution of the cold molecular gas, as described in the paper given in the references.
As can be seen in Table 1, all sources have been observed in either CO(1-0) or CO(2-1). For only two cases, data on more than one CO transition are available. The availability of more than one transition has been key to identifying the signatures of the interplay between the jet and the ISM, for example, by identifying gas with a high excitation temperature. The observations have spatial resolutions ranging from 0.2 to 2.0 arcsec, allowing us to spatially resolve the distribution of CO across the radio continuum emission in most cases.
In Section 3, we start by describing the new results from the NOEMA observations of 3C 305, 4C 52.37, and 4C 31.04, as well as providing a short description of the general properties of these three objects. In Section 4, we present a summary of the CO observations of the remaining objects in the sample, including the results obtained at other wavebands relevant to the final interpretation of the CO data.  (Murthy et al. in prep).

New
The data reduction followed the standard steps. The initial flagging of bad visibilities and the calibration were performed using the Grenoble Image and Line Data Analysis Software (GILDAS). Next, the two polarisations were averaged and the calibrated visibilities exported to uvfits format for further reduction in CASA. The data were initially self-calibrated on the continuum emission via a few cycles of imaging and phase-only calibration, followed by a round of amplitude-and-phase calibration. The line-free channels were then used to fit a firstorder polynomial to each calibrated visibility spectrum and subtracted from the uv data. Cubes were produced with a robust parameter of -1 and natural weighting in order to characterise the distribution and kinematics of the molecular gas, ensuring that sure low-surface brightness emission could be detected. The parameters of the final cubes are summarised in Table 2.   We de-redshifted the uv data to the systemic velocity of the target using the redshifts listed in Table 3 and obtained the final spectral cubes by imaging these datasets after averaging four channels. The velocity resolutions of the cubes are listed in Table 3.
The cubes were visually inspected for the presence of CO(1-0). However, to ensure that the analysis was systematic, they were also run through the Source Finding Application (SoFiA; [31]) 1 . This software performs spatial and velocity smoothing in order to optimise the search for emission.
In order to derive the molecular mass, we first estimated the CO(1-0) luminosity L CO(1−0) [32,33]: where z is the redshift, D L the luminosity distance in Mpc, and S CO(1−0) ∆v the integrated flux in Jy km s −1 of the CO(1-0) line. The molecular mass was estimated using where α CO is the CO-to-H 2 conversion factor in units of M (K km s −1 pc 2 ) −1 . The derived masses of the molecular gas are given in Table 3 using α CO = 3.6 M /(K km s −1 pc 2 ) [32]. This choice of α CO is based on the fact that, as discussed below, a large fraction of gas is regularly rotating. However, this is not the case for all the regions (in particular in 3C 305 and 4C 31.04, see below) and, therefore, a different (lower) conversion factor could be required. Thus, the values of the H 2 mass listed in Table 3 should be considered the upper limits. The 3 mm continuum images were obtained from the line-free channels and were made with the same weightings as those used for the line cubes.  The columns are (1) source name, (2) redshift, (3) molecular gas mass estimated, assuming α CO = 3.6 M /(K km s −1 pc 2 ). See text for details.
The width of the observing band allows us to explore the presence of the CN(1-0) lines. Interestingly, these were detected (in emission and absorption) in two of the targets (4C 52.37 and 4C 31.04). We did not detect CN(1-0) in 3C 305, possibly due to the limited sensitivity of these observations. The discussion of the CN lines is presented in Section 3.4.
In the following sections, we describe the results obtained from the CO(1-0) observations of the three targets.

3C 305
This radio galaxy and the multi-phase gas in which it is embedded have been the topic of several studies [30,[34][35][36]. The radio continuum morphology is peculiar, with two bright jets and hotspots, as well as lobes, which extend mostly perpendicular to the radio axis (see [30] and the references therein, and Figure 2). This structure suggests the presence of an interaction between the radio plasma and the ISM. Indeed, such an interaction has been traced, particularly on the E side, by the kinematics and properties of the ionised and H I gas, which show broad lines and outflows ( [35,36]). The physical conditions and kinematics of the warm molecular gas have been investigated by [34] using the Spitzer Infrared Spectrograph (IRS). This revealed multiple warm H 2 lines with FWHM above 500 km s −1 , indicating that the warm H 2 gas is also very turbulent in this source. In line with this and based on the X-ray morphology and spectral properties derived with Chandra, [30] argued that the extended X-ray emission was likely produced by shock-heating from interaction with the radio jet. Cold molecular gas was earlier detected in 3C 305 using single-dish observations with the IRAM-30m telescope [37,38]. Our new NOEMA observations have the angular resolution to trace the distribution of the gas. The molecular gas is distributed in a rotating structure roughly aligned with the radio axis and follows the warped dust lane, as shown in Figure 2 (orange contours). The total intensity (moment-0) image and the velocity field (moment-1) image overlaid with the contours from the 3 mm continuum emission are shown in the top panels of Figure 3. As can be seen in the moment-0 image (where the location of the radio core is also indicated), the total extent of the distribution of the molecular gas is about 17 arcsec, corresponding to 14 kpc. A brighter, inner region extending to ∼4 arcsec (about 3 kpc) can be observed. This is mostly co-spatial with the inner radio continuum emission, whereas the brightness decreases at the location of the hotspots (see below). The molecular gas extends on the W side to about 8 arcsec (6.6 kpc) following a warped structure. On the E side, most of the molecular gas extends to 4 arcsec, whereas more but fainter emission is detected up to 9 arcsec (7.4 kpc). The velocity field shows the overall rotation, which dominates the CO structure, and even the cloud east of 3C 305 follows these kinematics. However, at radii larger than ∼2 arcsec (1.6 kpc), in particular, on the E side, the gas distribution and kinematics of the molecular gas appear to be less regular. This is clearer in the position-velocity plots obtained along the major axis of the molecular gas along the radio axis (Figure 3 bottom). The dashed vertical lines indicate the locations of the hotspots in the radio continuum, showing the broad emission (>150 km s −1 ) seen, in particular, on the E side.  Top: 3C 305 moment-0 (left) and moment-1 (right) maps obtained using a mask to include only emission above the 4-σ level. This was made from the cube with robust weighting (see Table 3). The contours show the radio continuum map at a similar spatial resolution. Contour levels in both top panels are 0.2, 0.8, 2.8, and 10.0 mJy beam −1 . The location of the radio core is indicated. Bottom: position-velocity plot along the major axis (left) of the distribution of the molecular gas and (right) along the radio axis. The vertical dashed lines indicate the extent of the radio source. Contour levels in both position-velocity diagrams are 1, 2, 3, 4, 5, 6, 7, and 8 mJy beam −1 .
We derived a total molecular gas mass of 2.3 × 10 9 M using a Galactic CO-H 2 conversion factor of α CO = 3.6 M /(K km s −1 pc 2 ). This is consistent with that derived by [38], who found 2.1 × 10 9 M of molecular gas. However, since we also found that in some regions, the gas showed high-velocity dispersion (pointing to optically thin gas, see below), a lower conversion factor may be more appropriate for those regions. In this case, the estimated molecular gas mass represents the upper limit.
Interestingly, the location of the broad emission is also where H I and ionised gas outflows were detected and could represent the molecular counterpart of these outflows. The broad CO emission is too faint for a proper characterisation of the kinematics of this gas. Nevertheless, if we assume that this gas represents the counterpart of the H I outflow, some useful constraints on the properties of a possible molecular outflow can be derived. We estimated the upper limit of the mass of this potential molecular outflow by considering a 3-σ signal over 400 km s −1 (i.e., a width comparable to the outflow observed in H I). Using a conversion factor of α CO = 0.34 M /(K km s −1 pc 2 ), typical of the optically thin regime found in outflows, we derived an upper limit of M out H 2 ∼1.5 × 10 7 M . The corresponding mass outflow rate wasṀ∼4 M yr −1 . Considering the large uncertainties due to the assumptions made, this is consistent with the H I mass outflow rate estimated in [36], suggesting that the molecular outflowing component, if present, is not large.
A direct comparison between the H I and CO outflows is complicated by the fact that the H I outflow has been observed in the absorption and thus can only be traced against the regions where the radio continuum is detected. The molecular gas seen in the emission suggests that the situation may be more complex. For example, we note that the faint emission of CO gas with broad profiles is mostly located east of the peak of the continuum (indicated by the dashed line in Figure 3). Additionally, it is interesting that the total intensity image (Figure 3 top left) suggests that the molecular gas (especially on the E side) avoids the brighter part of the radio lobe. Thus, in the case of 3C 305, we may be seeing the radio lobe in the process of clearing the molecular gas by pushing it aside. This could be similar to what was seen in PKS 0023-26. There, the molecular gas with broad, albeit not extremely broad, profiles was seen wrapping around the lobe (see [39] for details). We revisit this in the discussion in Section 5.

4C 52.37
4C 52.37 is a radio galaxy that is included in the CORALZ sample of compact and young radio sources selected by [19]. In the VLBI sub-arcsecond radio morphology, it has been classified as a compact symmetric object [40].
The optical spectrum of 4C 52.37 is characterised by strong, broad permitted, and forbidden emission lines. This radio galaxy is known to be embedded in an H Irich medium [24,41], where the H I appears kinematically disturbed. Broad, blueshifted H I absorption was found covering more than 600 km s −1 [24]. The broad H I absorption is mostly concentrated in the central regions (i.e., the central ∼100 pc) of the radio source [41]. The morphology of the radio continuum on VLBI scales is shown in the inset in Figure 4, whereas the mm continuum from the NOEMA observations is unresolved (contours in Figure 4).
For this object, our CO(1-0) observations are of relatively low resolution with respect to the size of the radio source. This is because no previous observations of the molecular gas were available and, therefore, a conservative approach was used for our detection experiment. Higher resolution observations are needed for a more complete assessment of any jet-molecular gas interaction.

kpc
East West Figure 4. Left/Centre: 4C 52.37 moment-0 and moment-1 maps obtained using a mask to include only emission above 4-σ level. This was made from the cube with robust weighting (see Table 3). Contours show the radio continuum map at a similar spatial resolution; contour levels are 2.5, 5.0, 7.5, and 10 mJy beam In the NOEMA observations, we detected CO(1-0) located in a regular disc of 12 × 10 9 M (assuming the Galactic value α CO = 3.6 M /(K km s −1 pc 2 )) and distributed over a total extent of 10 arcsec (∼19.4 kpc). The gas reaches larger radii (∼6 arcsec, 11.7 kpc) on the SE side. The moment maps and the position-velocity plot are shown in Figure 4. In particular, the position-velocity plot is characteristic of a regularly rotating disc. The molecular disc is oriented in PA 129 • , which is a 40 • offset (in projection) from the radio axis ( [41]; see the inset in Figure 4). Considering the axis ratio of the integrated CO emission, the disc has a fairly high inclination (i.e about 60-70 • ) but it is not edge-on.
As mentioned above, an H I outflow was detected in absorption against the central regions (∼100 pc) of this source. Using a similar approach as for 3C 305, we can estimate a limit to the molecular outflow that could have been detected in our observations. We obtain about 8 × 10 7 M for a 3-σ profile covering an approximately 400 km s −1 width. If the molecular outflow were located where the H I outflow is found (between 60 and 150 pc in the E radio lobe), the mass outflow rate of the molecular outflow would be less than ∼30 M yr −1 . Compared to what was detected in H I (∼4 M yr −1 ; [41]), this suggests that a massive molecular outflow is not present but our observations are not sensitive enough to trace a more modest outflow.

4C 31.04
This is a well-studied radio galaxy with a wealth of ancillary data [28,42,43]. The VLBI radio continuum emission shows a core and two lobes separated by only ∼100 pc (inset in Figure 5 middle). The lobes are very asymmetric, with the western lobe showing diffuse emission unlike the eastern lobe, which has a bright hotspot. The host galaxy was observed with HST by [42] and shows an edge-on circum-nuclear dust lane with a prominent warp ( Figure 5 left). We observed 4C 31.04 both in H I using VLBI and in CO(1-0) with NOEMA. The full presentation and discussion of the data can be found in the paper by Murthy et al. (currently being prepared), whereas here, we provide the highlights of the CO observations. Molecular gas single-dish observations were carried out in CO(1-0) and CO(2-1) by [38], who detected 6.16 × 10 8 M of H 2 in emission and a strong absorption line in CO(1-0). Using IRAM Plateau de Bure (PdB), [45] detected a circumnuclear disc of about 1.4 kpc in size in HCO+, both in emission and absorption. Although the kinematics of the gas are dominated by rotation, these data show evidence of distortions and non-circular motions, which suggest the disk is not in a dynamically relaxed state and that it could provide the reservoir for fuelling the radio source. The presence of a jet-ISM interaction is further supported by the observations of warm molecular gas [43]. The detection of [FeII] lines in the inner region has been described by these authors as tracing shocked gas ejected from the disc plane by a jet-blown bubble (300-400 pc in diameter), whereas the warm H 2 emission traces shock-excited the molecular gas in the interior ∼1 kpc of the circumnuclear disc [43].
Our NOEMA CO(1-0) observations revealed a disc structure with a full extent of about 9 arcsec (10.4 kpc) and a mass of 2.4 × 10 9 M (assuming the Galactic value α CO = 3.6 M /(K km s −1 pc 2 )). This structure is shown in the middle in Figure 5. This disc follows the outer region of the dust lane seen by HST, as shown by the similar position angle. This complements the 1.4 kpc-sized HCO+ disc presented by [45], which instead traces the inner regions. The difference in position angle between this inner disc and the one observed in CO(1-0)confirms that the warped structure observed in the dust lane is also traced by the cold molecular gas.
Strong absorption was seen against the 3 mm continuum, which is unresolved at our resolution. The absorption has a full width at zero intensity (FWZI) of ∼250 km s −1 and is blueshifted with respect to the systemic velocity, as can be seen on the right in Figure 5. Interestingly, similar properties were also found for the H I absorption observed with VLBI (Murthy et al., in prep), which shows an FWZI of ∼ 300 km s −1 . However, although the CO continuum emission is unresolved (the size of the radio source is much smaller than the resolution of the present CO(1-0) observations), the VLBI resolved the continuum emission and allowed us to locate the absorption and derive the resolved kinematics of the H I gas. The H I absorption is clearly extended, covering most of the radio lobes, and blueshifted compared to the systemic velocity of the host galaxy. This indicates that it must be the result of kinematically disturbed gas not following the regular rotation of the large-scale disc. As mentioned above, the study of the warm molecular gas also suggests the presence of a strong interaction and, in particular, an expanding gas shell driven by the radio jets interacting with the ambient ISM. Thus, both the H I and the cold molecular gas appear to be part of this structure. In addition, for this target as for 4C 52.37, higher-spatial-resolution CO observations are needed to follow the molecular gas around the radio source in more detail.

CN(1-0) Lines
The broad band of NOEMA allowed us to investigate the presence of the cyanide radical (CN) lines simultaneously with the CO(1-0) lines. We detected the CN lines in two of the three objects, 4C 52.37 and 4C 31.04, and we show the integrated spectra in Figure 6. In the two sources, the CN lines exhibit a mix of absorption and emission. In 4C 31.04, the CN feature has the same structure as the CO, both in emission and absorption. However, the CN(0-1) absorption is deeper than that of the CO because the CN molecule has a much larger dipole moment than CO, making the critical density of the corresponding CN(1-0) line larger by a factor of 10-100 than that of the CO(1-0) line. In both sources, two sets of lines are seen in CN. In 4C 31.04, the two absorption lines have a peak line-strength ratio of about 2:1. The deepest absorption line is due to the CN (N = 0-1, J = 3/2-1/2) group of hyperfine structure lines (unresolved here at a resolution of 24 km/s, centred at a rest frequency of 113.49 GHz), and the weaker feature at a velocity of +750 km/s is the CN (N = 0-1, J = 1/2-1/2) group (centred around 113.16 GHz). This feature was also detected, although weaker, in 4C 52.37.
CN molecules are primarily produced by photo-dissociation reactions of HCN in the moderately dense (≈10 3 cm −3 ) molecular phase of the ISM [46]. Thus, CN emission lines mostly trace molecular gas irradiated by strong UV radiation fields. We estimated the integrated CO/CN ratios for the emission and absorption features using the brighter CN line, with the caveat that the absorption could be partly filled by emission. In 4C 52.37, we found a ratio of CO(1-0)/CN = 4 for the total integrated profile. In 4C 31.04, CO(1-0)/CN = 3 for emission and CO(1-0)/CN = 0.8 for absorption. A quantitative exploitation of those line ratios can be found in an as-yet unpublished paper but here, we qualitatively discuss why those ratios are significantly lower than those in nearby galaxies, which have, on average, CO/CN ≈ 10, with internal variations of a factor of 3-10 depending on the environment [47,48].
Models suggest that the CN/CO abundance and line flux ratios are enhanced in lowdensity, irradiated molecular shocks (such as the CH+/CO ratio, [49,50]); high-velocity dispersion diffuse molecular gas [51]; or X-ray-dominated regions [52]. For instance, high column densities of CN have been found close to AGN [47,53] and elevated CN/CO line flux ratios have been detected in the outflow of Mrk 231 [54] or from an absorbing merging system tidal tail towards the background quasar G0248 [55]. Given the conditions originating from the kinematically disturbed gas in our sources, it is likely that one or a combination of some of these processes could be responsible for the elevated CN flux relative to CO(1-0) in 4C 52.37 and 4C 31.04. A more detailed morpho-kinematical study of line ratios compared with chemical and dynamical modelling could identify the dominant mechanisms.
In 3C 305, we did not detect the CN lines. Considering the rms noise in the cube of ∼0.5 mJy beam −1 and the CO signal between 5 mJy (NE) and 4 mJy (SW), this gives a 3-σ lower limit for the ratio of CO/CN of ∼3. Thus, this does not allow us to draw strong conclusions about the conditions in this object.

Results from the Other Young Radio Galaxies in the Sample
Here, we briefly summarise the main results obtained from the CO observations of the other objects in the sample, before combining the results of all sources in Section 5 and outlining the picture obtained so far. For a detailed discussion of the results, we refer the reader to the published papers listed in Table 1.
The first object we discuss is B2 0258+35 (L 1.4 GHz = 2.1 × 10 23 W Hz −1 ), highlighting the impact of relatively low-power radio jets, as presented in [16]. NOEMA CO(1-0) observations have shown that the circumnuclear molecular gas is highly turbulent and a fast (FWHM∼350 km s −1 ) jet-driven outflow of ∼ 2.6 × 10 6 M has been detected (Figure 7 reproduced from [16]). This outflow dominates the kinematics of the molecular gas in the central kpc and makes up ∼75% of the total gas in the nuclear region. At this rate, the jet will deplete the kpc-scale molecular gas reservoir on a relatively short time scale of ∼2 × 10 6 yr. This finding confirms the predictions of the numerical simulations described in Section 1 and shown in Figure 1. The results on B2 0258+35 confirm that low-power jets (in the case of this source, of about 10 44 erg s −1 [16]) have enough energy to drive outflows. Interestingly, this was also seen in the Chandra X-ray observations of this source, which revealed extended X-ray emission both along and orthogonal to the jet [56]. At the end of the main kpc jet, lower-energy X-ray emission is seen coinciding with the region where the molecular gas shows high turbulence and fast outflow motions. This suggests that the hot ISM may be compressed by the jet and that the molecular outflow can result from more efficient cooling, which is in agreement with the predictions of the numerical simulation of the jet-ISM interaction for low-power jets.
It is interesting to note that for one of the best examples of a low-power radio jet impacting the molecular gas, IC 5063 [57], the Chandra X-ray observations confirm the presence of strong shocks through the detection of localised Fe XXV emission in the region where the jet-ISM interaction occurs (see [58] for details). Thus, these very different phases of the gas (ranging from cold molecular to hot gas) appear to provide a consistent picture of the interplay between jets and the ISM.
We also expanded our study of the molecular gas to young radio galaxies with highpower jets; the most extreme cases are PKS 1549-79 [59] and PKS 0023-26 [39]. In addition to powerful radio sources, these two objects also host a quasar-like optical AGN. Emitting both radiative and mechanical energy, the impact of the AGN is expected to be quite prominent in these sources.
In PKS 1549-79, using ALMA, we found the most massive outflow of molecular gas in the sample (∼650 M yr −1 ; [59]). However, the outflow is limited to the inner 200 pc of the galaxy despite the presence of a powerful jet (about 300 pc in size) in combination with a powerful quasar AGN, with bolometric luminosity between 6 and 9.6 × 10 45 erg s −1 [21]. A circumnuclear disc of M H 2 = 2.6 × 10 8 M was also observed and appears to co-exist with the outflow. This would suggest a depletion time of ∼10 5 yr, although a tidal tail of molecular gas, which was still in the process of falling in, was observed and may help to refuel the disc (and likely the AGN; see [59] for details).
PKS 0023-26 is also radio powerful but more extended, i.e., it is one of the largest radio galaxies in the sample [39]. The radio emission of this radio galaxy is a few kiloparsec in size, thus probing interesting intermediate scales where the young jet may start to break free from the inner dense regions of the ISM. The presence of direct interactions between the jet and the ISM was seen in the broad profiles of the CO observed in the inner kpc region. On larger scales, the molecular gas instead appears to avoid the radio lobes and wrap around them. Although the kinematics of the gas are disturbed, the profiles are less broad than in the centre and the gas appears to be affected not by the direct interaction with the jet but by the expanding cocoon surrounding the radio source, likely dispersing and heating pre-existing molecular clouds (see [39] for details).
These results suggest an evolutionary sequence, similar to that predicted by the numerical simulations (e.g., [5]) in which the impact of the jet on the surrounding medium changes with its expansion. We revisit this in Section 6.

The Picture So Far
From studying the cold molecular gas around eight young radio galaxies, we found a variety of physical and kinematic properties of the observed CO structures, suggesting a complex interplay between jets and cold molecular gas. This was not unexpected given that it was shown by numerical simulations (e.g., [7]) and several parameters play a role in defining the impact of such an interaction. We also found some interesting trends, which can be compared with the predictions of the simulations.
In most objects studied (with the possible exception of PKS 0023-26 and the central region of B2 0258+35), we saw a large fraction of the molecular gas distributed in rotating structures. The full sizes of these structures cover a broad range from sub-kpc to about 20 kpc. Similar sizes of molecular discs were also seen in evolved, large-scale radio galaxies and radio-quiet early-type galaxies (e.g., [60][61][62][63] and Tadhunter et al., currently being prepared). Assuming it would take at least a few rotations to form such structures, the discs we observed must have formed between a few ×10 6 yr ago for the discs of hundreds of pc and ∼10 8 yr ago for the largest (a few kpc) discs. They were, therefore, in place before the radio source (re-)started. This suggests that, although these structures can be impacted in various ways by a young jet (see below), they are seldom completely disrupted by the onset and evolution of the jet, as illustrated by the simulations in Figure 1.
Our observations show, as seen previously for other phases of the gas, that young ( < ∼ 10 6 yr) radio jets are able to drive fast and massive gas outflows (see IC 5063, B2 0258+35 and PKS 1549-79; [16,57,59] for the clearest cases). However, the results from the cold molecular gas show, more clearly than in the case of the ionised gas, that low-power jets (i.e., P jet < 10 45 erg s −1 ) are able to drive such massive outflows. This is the case for IC 5063 and B2 0258+35. This confirms the predictions of the simulations, which show that, although low-power jets can temporarily remain trapped in the ISM, they can accumulate enough energy during this time to break through the gas and drive outflows [7]. Because they are less powerful, the entire process takes longer than for powerful jets but can be extremely effective. This is relevant because it highlights the relevance of the feedback of low-power sources, which are more common than powerful sources.
However, all the observed molecular outflows are limited to the central (at most kpc) regions. This is also the case for objects (such as the two most powerful sources in the sample, PKS 1549-79 and PKS 0023-26) where both a strong radio and an optical AGN are present and, therefore, the AGN impact is expected to be high.
In two cases, the depletion time of the outflows could also be derived. For B2 0258+35, the derived mass outflow rate (between 5 and 10 M yr −1 ) implies that the kiloparsec-scale molecular gas reservoir will be depleted in ∼2 × 10 6 yr (see [16] for details). Even shorter is the depletion time found for PKS1 549-79, where the massive outflow would deplete the central disc on a time scale of only ∼10 5 yr. Although we are limited by the small number of objects studied, these results suggest that the conditions for developing an outflow could last only a relatively short time, perhaps explaining the fact that molecular gas outflows are not observed ubiquitously. The outflow phase can, however, be recurrent in the cycles of activity that are known to occur during the life cycles of radio galaxies (see [64] and the references therein).
In addition to outflows, we found disturbed and highly turbulent gas in the regions co-spatial with the radio emission in most of the objects. This is another sign of the effect of the expansion of the young radio lobes into the surrounding medium.
A further indication of the impact of the expanding jet can be found in the physical conditions of the gas. These appear to be different in the kinematically disturbed gas from those in regularly rotating gas. We observed this, for example, in the CO line ratios in IC 5063 and PKS 1549-79 (the two objects where multiple transitions are available; see [59,65,66]). The line ratios from different CO transitions were found to be substantially different in the kinematically disturbed gas from what is typically found in relaxed galaxy discs, with the higher transitions being relatively much brighter in the disturbed gas. The observed line ratios imply that the molecular gas near the AGN is optically thin and has elevated excitation temperatures (T ex > 50 K) as a result of the impact of strong shocks. Similar results have been found by [67] in the case of the young radio source hosted by NGC 3100. Using multiple CO transitions, these authors found a high-excitation molecular gas component in the region close to the radio jet, whereas the gas becomes progressively less excited at increasing distances from the radio structure. The fact that in NGC 3100 only mild signatures of the interaction can be seen in the kinematics [67] suggests that the line ratios can be a more sensitive indicator of the presence of jet-ISM interactions.
Thus, to obtain a complete picture of the impact of the radio plasma on the ISM, not only are the kinematics of the gas needed but also observations of multiple transitions (and/or multiple molecules) in order to derive the conditions of the gas such as the excitation temperature. This is also shown by the detection of CN lines in two of our targets. The broad band of NOEMA allowed us to not only detect these lines but also measure the CN/CO fraction. The high fraction detected compared to that typically found in quiescent gas discs can be explained as the result of the impact of the AGN on the conditions of the gas. More detailed investigations are required in order to fully quantify the implications.
One of the goals of this project was to investigate whether the interaction between the radio jet and the ISM changes as the jet expands. Our results suggest that this change indeed occurs. For the largest objects in our sample (PKS 0023 -26 and 3C 305, both a few kpc in size), the molecular gas appears to avoid the brighter part of the lobes and instead wraps around them. This is particularly clear in PKS 0023-26 in the channel maps shown in Figure 5 of [39] and a similar trend was seen in 3C 305.
Because of this change (albeit based on only a few cases), we propose a possible evolutionary scenario (see also Figure 11 in [39]). The fast and massive outflows are limited to the very inner regions close to the centre, where the newly emerged jet interacts directly with the gas clouds while opening its way in the surrounding ISM. When the jets reach the kpc radii, this is replaced by a milder, more gentle interaction driven by the expanding (both along and orthogonal to the jet) cocoon created as a consequence of the initial direct jet-cloud interaction. This results in the dispersing and heating of the ISM. These findings, which are in line with predictions from numerical simulations of jets interacting with a clumpy medium, also suggest that the effect of expanding shells of shocked gas-reminiscent of the "maintaining mode" of feedback associated with X-ray cavities in cluster galaxies-could provide an important way for radio plasma jets to affect the galactic ISM. In our study, we have seen several signatures of interaction between the radio plasma jets and the ISM in the vast majority of young radio sources located in a gas-rich environment. This kind of interaction may, in turn, also affect the general evolution of the radio jets. A relationship between the turnover frequency of the radio spectrum and the size of the radio source has already been suggested in early studies of young radio galaxies (i.e., [17,68]). This relationship has been explained by the turnover frequency being due to synchrotron self-absorption or freefree absorption due to a dense medium and it evolves as the source grows. It is interesting to see that the turnover frequencies of the sources studied here (in some cases, derived by adding new low-frequency points, e.g., from GLEAM or LOFAR, see, e.g., [23,39]) are all located below the correlation-size turnover frequency from the literature, as shown in Figure 8. The most likely explanation for this is that radio sources that strongly interact with a rich medium are affected in the speed of their expansion, resulting in a smaller size given the turnover frequency observed. This would suggest that the expansion of the radio lobe has been slowed down due to the interaction with the rich medium.
As a final remark, it is interesting to note that the only object where we observed clear signs of molecular gas in the process of falling into the SMBH was the smaller radio source in the sample, PKS 1718-63. With only one object available, we can only speculate that the very young age of the jet means that the jet-ISM interaction is yet to develop. Thus, in this object, we can more clearly see the gas that is instead in the process of feeding the AGN before the disentangling of the two processes (feeding and feedback) becomes too complex.

Discussion and Future Work
We explored the jet-ISM interplay for a small sample of young radio galaxies in which the newly born jet is expanding into the rich ISM of the central regions of the host galaxy. We found that in the vast majority of the objects (with perhaps one exception), the gas is kinematically disturbed in the region co-spatial with the radio galaxy, showing the coupling between the jets and the ISM.
In only a small fraction of the sources were the observed kinematics of the molecular gas extreme. However, the finding that radio jets and, in particular, low-power jets, are able to produce fast and massive outflows, makes this type of AGN relevant for cosmological simulations, even on galactic scales. Indeed, some recent cosmological simulations have started including more realistic descriptions of the impact of the jets [9] despite the huge challenges involved. We further observed that gas outflows are limited to the very inner regions of the galaxy. Thus, the outflows observed (even in powerful AGN) do not seem to be large or extended enough to provide feedback and affect star formation on galactic scales.
Interestingly, in the largest radio sources (a few kpc in size), we observed the gas avoiding the region with the radio continuum emission and wrapping around the radio lobes. Thus, we suggest an evolution in the type of interaction between radio plasma and gas. Outflows from the direct jet/cloud interaction in the clumpy medium are occurring in the inner region (few hundred parsecs/sub-kpc), whereas on larger scales (larger than a kpc), the impact is mediated by the cocoon of the shocked gas created by the jet-cloud interaction; this cocoon pushes the gas aside and heats it. Consistent with this is that a number of our targets observed with Chandra show signatures of X-ray gas thermally heated by the interaction between the radio plasma and the ISM.
If this scenario is correct, we argue that the impact of the radio plasma will continue on galactic scales as the jets expand. The impact of galactic-scale jets (GSJ; [69,70]) has already been investigated in a handful of objects, with observations looking at the distribution of X-ray gas. Using the Chandra observations of these sources, [11][12][13][14] have found the presence of structures of hot gas thermally heated by shocks from the jets, as shown by the morphological correspondence and distribution around the radio lobes. Thus, in this phase, the role of the radio plasma could mostly be preventing the cooling of the gas, similar to what occurs in the hot medium in clusters by the creation of X-ray cavities. A detailed study of the impact of GSJ on both the hot and cold/warm gas could be the way forward to obtain a full overview of the role of the radio plasma jets in feedback.
The trends we have found for a small number of objects need to be confirmed by observations of larger samples. Although this is a demanding task requiring high spatial resolution and sensitive observations with oversubscribed telescopes such as ALMA and NOEMA, the continuously expanding ALMA archive will provide a growing database, which will hopefully soon allow for the expansion of the present results. Funding: This research received no external funding.

Data Availability Statement:
The data from the NOEMA archive and the final calibrated data, images, and cubes are available on request to the authors.