Numerical simulations of colliding jets in an external wind: application

The radio galaxy 3C 75 is remarkable because it contains a pair of radio-loud active galaxies, each of which produces a two-sided jet, with the jet beams appearing to collide and merge to the west of the galaxies. Motivated by 3C 75, we have conducted three-dimensional hydrodynamic simulations of jet collisions. We have extended previous studies by modelling the physical properties of the cluster atmosphere, including an external wind, and using realistic jet powers obtained from observational data. We are able to produce a morphology similar to that of 3C 75. The simulations imply that direct contact between the bulk jet ﬂows on the west of the source is required to produce a morphology consistent with 3C 75. We quantify how the merging jets decelerate, how the wind deﬂects the jets and cocoons, the entrainment of intra-cluster material into the cocoons, the cocoon energetics, and how the jet interactions generate enstrophy. By comparing simulations of pairs of two-sided jets with those of single two-sided sources, we determine how the interaction between two bipolar jets changes their evolution. The unprecedented sensitivity and angular resolution of upcoming observatories will lead to the detection of many more complex sources at high redshift, where interacting jets are expected to be more numerous. The morphology of these complex sources can provide signiﬁcant insight into the conditions in their environments.


I N T RO D U C T I O N
The 3C 75 radio source is hosted by a pair of merging galaxy clusters that contain two prominent radio-loud active galactic nuclei (AGN). Remarkably, the jets from these AGN interact with each other and appear to collide (see Fig. 1). The two active nuclei that make up the 3C 75 radio source are located in galaxy NGC 1128, in the Abell 400 galaxy cluster. The redshift of 3C 75 is z = 0.023 (Owen et al. 1985;Hudson et al. 2006), which corresponds to a luminosity distance of 101 Mpc and an angular scale of 1 = 0.47 kpc. The nuclei of 3C 75 are offset from the X-ray emission-weighted centre of the cluster by 1.68 arcmin (≈48.9 kpc; Hudson et al. 2006). The supermassive black holes of the AGN are separated by a deprojected distance of 8 kpc (7 kpc in projection when observed with an inclination angle of 20 • for the orbital plane of the two nuclei; Yokosawa & Inoue 1985;Hudson et al. 2006). Jet collisions are rare events and 3C 75 presents a rather striking example as jets on both E-mail: g.musoke@uva.nl sides of the source display significant interaction. A particularly strong interaction occurs ∼10 kpc north-west of the nuclei, where the two jets on this side of the source collide, appear to merge, and then continue as the northern plume. The eastern side of the source comprises a comparatively weaker encounter where the jet plumes appear to interact with each other and display significant curvature. The large-scale morphology of the source is reminiscent of Wide-Angle-Tailed (WAT) sources, displaying plumes that are preferentially bent towards the north-east of the source. The jets are bent due to the relative motion between the cluster gas and the source. Due to the non-thermal nature of relativistic jet spectra, fundamental jet properties are difficult to discern. However, the interactions of jets with their environments can provide an insight into these jet properties.
The first well-resolved radio observations of 3C 75 were conducted by Owen et al. (1985) in the VLA Abell cluster snapshot survey. Observations of the source prior to this (e.g. Fomalont 1971) were not well resolved and did not show the true nature of the source due to its low declination (Owen et al. 1985). Early Xray observations of 3C 75 presented by Forman & Jones (1982) Figure 1. Adaptively smoothed, exposure map corrected 0.5-7 keV Chandra X-ray image (orange) of 3C 75 with contours of 4.8 GHz VLA radio map overlaid (blue). The radio source contains four jets that interact with each other on both sides of the source. A particularly strong interaction is observed on the north-west side of the source, where the jets appear to collide and merge, though the true nature of the interaction remains unknown. The jets could be colliding head on, or glancing off each other.
were conducted with the Einstein X-ray observatory. These X-ray observations revealed the elongation of the cluster gas along the northeast-southwest direction (Beers et al. 1992), the line towards which the jets are bending. Consequently, Beers et al. (1992) suggested that the alignment between the jet bending and cluster elongation implied that the jet bending observed in 3C 75 was due to the relative motion of the source as it moves through the intracluster medium (ICM). A number of authors (e.g. Begelman, Rees & Blandford 1979;Balsara & Norman 1992;Canto & Raga 1995;Morsony et al. 2013;Gan et al. 2017;Jones et al. 2017) have explored the deflection of single two-sided AGN jets by an external wind. 1 The general consensus is that AGN jets bend due to the ram pressure of a cross wind as a result of the motion of the source relative to the surrounding ICM. The plausibility of this scenario for 3C 75 was confirmed through the kinematic analysis of the host galaxy cluster Abell 400 by Beers et al. (1992), who concluded that the galaxy cluster is likely currently undergoing a merger of two primary subclusters of similar mass. In this picture, the collisional velocity of the merging subclusters must be large in order to cause the large-scale bending of the radio jets.
Using the Chandra X-ray observations alongside 4.5-GHz and 330-MHz VLA radio data, Hudson et al. (2006) concluded that the AGN are contained in a low-entropy cluster core and are a bound system from a previous merger. Though the radial velocities of the two nuclei differ by more than 400 km s −1 and their angular positions differ by ∼16 arcsec, the nuclei are likely bound. This is suggested by the high mass of the combined galaxy and the fact that the jets ejected by each AGN are observed to bend in a similar direction -the latter suggests that both AGN must have the same relative motion with respect to the ICM (Hudson et al. 2006). The host galaxy NGC 1128 is believed to lie in one of two dynamically distinct galaxy subsystems within the cluster Abell 400. Hudson et al. (2006) also conclude that the Abell 400 cluster is currently undergoing a merger, in agreement with the dynamical analysis of Beers et al. (1992). This cluster merger scenario is further supported by the absence of a cool core in Abell 400 (Jones & Forman 1984;Hudson et al. 2006), which is consistent with a merging system, and evidence of a shock south of the cluster core (Hudson et al. 2006).
The X-ray observations presented by Hudson et al. (2006) also revealed the presence of an elliptical galaxy, referred to as A400-41 (PGC 011193), that is located ∼40kpc to the east of the nuclei. A400-41 is almost at the same redshift (within 400 km s −1 ) as NGC 1128. The authors found that A400-41 is bright in X-rays but displays no radio emission. Consequently, Hudson et al. (2006) assume that the emission associated with A400-41 is thermal in nature and thus the elliptical galaxy is not thought to be active. A400-41 appears to interact with 3C 75, with the radio emission of the eastern jet from the southern AGN appearing to bend suddenly as if deflected by the galaxy. Simulations conducted by Wang, Wiita & Hooda (2000) have shown that it is possible for jets in some circumstances to be deflected by galaxies. By estimating the kinetic luminosity of the eastern jet of the southern AGN that interacts with A400-41, and estimating the minimum jet flow speed required to ensure that the pressure exerted by the jet on the galaxy is less than the average pressure of the galaxy, Hudson et al. (2006) conclude that it is possible for A400-41 to cause the deflection if the jet is light and supersonic.
While the large-scale bending of the source can be attributed to the relative motion between the nuclei and the ICM, it is not clear what causes the morphology of the interacting jets on the northwest side of the source. This jet-jet interaction could be interpreted as (1) the jet trajectories crossing, with the jets twisting around each other in large-scale helical motions, 2 (2) the jets bending in tandem with each other, or (3) the jets colliding and merging together. Determining which of these interpretations, if any, is most appropriate for 3C 75 through observational data alone is difficult. It is unclear whether these north-western jets actually collide, merely glance off each other, or simply do not come into contact at all. While it is possible that the jets propagate one behind the other and appear only to collide or interact in projection, we note that the radio emission at the location of the putative interaction is brighter than expected from projection alone. This supports the interpretation that the north-western jets are physically interacting with one another.
A number of models have previously been suggested to explain the large-scale morphology of 3C 75 and particularly the intertwining structure of the jets on the north-western side of the source. Yokosawa & Inoue (1985) conducted early numerical simulations of intertwining jets, focusing on recreating the morphology observed in 3C 75. They concluded that the wiggling of the north-western jets could be explained by their proposed 'binary-orbiting radio jet' model. In this model, the radio lobes associated with each of the jets interact with each other, while their cores orbit each other with an orbital period of ∼10 8 yr and an orbital radius of 8 kpc. The model assumes that the density of the jets is much higher than the density of the surrounding ICM. The jets are assumed to consist of independent plasmoids that are ejected from their respective nucleus, propagate ballistically, and are decelerated by the ram pressure of the surrounding ICM. Yokosawa & Inoue (1985) suggest that the large-scale radio morphology of 3C 75, in which the ends of the diffuse plumes appear to be swept back in a direction nearly perpendicular to the direction of the ICM wind, could be attributed to the effects of buoyancy in the gravitational field of the cluster, which may dominate over ram pressure at larger distances from the cluster centre (Sakelliou, Merrifield & McHardy 1996).
Lastly, from their simulations, Yokosawa & Inoue (1985) estimate that the source must move through the ICM with a relative speed of 1120 km s −1 to produce the observed large-scale jet/plume bending. However, the binary-orbiting radio jet model proposed by Yokosawa & Inoue (1985) does not explain the absence of helical features, or the spiralling of the jets, on the south-eastern side of the source. If the helical topology of the jets on the north-western side of the source is due to the binary motion of the nuclei, then the orbital motion should be dynamically important and large-scale helical features should be observed on both sides of the source. The latter is not the case for 3C 75.
More recently, Hudson et al. (2006) presented a 22-ks X-ray observation of the source using the Chandra Advanced CCD Imaging Spectrometer (ACIS) instrument. The observation was conducted as part of the Highest X-ray Flux Galaxy Cluster Sample (HIFLUGCS) follow-up program (Reiprich & Böhringer 2002). Hudson et al. (2006) estimated a relative motion between the source and the ICM of ∼1200 km s −1 , with a Mach number M ∼ 1.4 (given a cluster sound speed of 850 km s −1 ), by noting the possible presence of a shock to the south-west of the cluster core. In agreement with Yokosawa & Inoue (1985) and Beers et al. (1992), Hudson et al. (2006) also noted that the intracluster 'wind', given by the relative motion of the source as it moves through the ICM, is likely the cause of the initial bending of the jets and that this wind, alongside turbulence, is also likely the cause of the jet/plumes bending on scales of 10's of kpc.
Another possible explanation for the large-scale radio morphology of 3C 75 is given by Achterberg (1988), who showed that the morphology can be caused by the direct interactions of jets in close proximity to each other. Achterberg (1988) expounded two models, one based on the electromagnetic interactions of the jets and the other governed by hydrodynamic interactions. In both models, the jet-on-jet interactions were considered in the fire hose limit, the underlying assumptions of which are that the cross section of the jet is not significantly distorted, the cross section varies slowly only with distance along the jet, and the large-scale motions of the jet are well represented by the motion of the jet axis such that the global movements of the jet correspond to that of a tensile string. In the electromagnetic model proposed by Achterberg (1988), the complex radio morphology can be produced by currentcarrying jets interacting electromagnetically in the presence of an external magnetic field. In this model, the jets are embedded in a common cocoon and 'feel' each other's magnetic fields, causing their trajectories to bend and resulting in a helical topology sustained via the Lorentz force. In the hydromagnetic model, one jet traverses the wake of the other produced by the relative ICM flow with a component perpendicular to the axis. The wake induces wiggles in the trajectory of the jet that moves through it due to the action of fluctuating drag and lift forces that result from fluctuations in the vorticity and velocity of the surrounding medium. The large-scale bending of the source is due to a high Reynolds number cross-flow of the ICM, which results from the relative motion of the host galaxy moving through the ICM.
Both the hydrodynamic and electromagnetic models proposed by Achterberg (1988) require that the jets are underdense with respect to the density of the surrounding ICM. For the electromagnetic case, Achterberg (1988) concludes that ratio of the jet density to the ambient density must be less than or equal to 0.01, in order for instabilities to act slowly enough that the helical structure on the north-western side of the source survives throughout the lifetime of the source. For the hydrodynamic case, this ratio is required to be less than or equal to 0.05 so that the lift and drag forces are strong enough to cause the wiggles observed in the jets.
In Molnar et al. (2017), we performed three-dimensional (3D) hydrodynamic simulations of colliding jets, motivated by the morphology of the north-western side of 3C 75. The simulations were performed with various jet velocities and spatial separation, or collisional impact parameter, between the jets. We define the collisional impact parameter as the offset in the z direction between two jets aligned in the x, y plane. Throughout the following discussions, the impact parameter between the colliding jets is given in units of the jet nozzle diameter and is measured from the nozzle centres (see Sections 2 and 5.3 for more details). In Molnar et al. (2017), we aimed to explore how interactions between binary twosided jets are influenced by the jet velocity and collisional impact parameter in order to determine which scenarios may favour the morphology observed in 3C 75.
In Molnar et al. (2017), we found that the behaviour of binary twosided jets following an interaction (or collision) between them can be grouped into two categories -bouncing and merging jets. In the bouncing jet category, the jets bounce or glance off each other and retain their individual identities following the interaction/collision. In the merging category, only one jet emerges from the jet collision. We found that jets with similar velocities and non-zero impact parameters tend to result in bouncing behaviour, while large velocity contrasts between the colliding jets and small impact parameters favoured merging jets. In addition, we found that in some cases the collision can enhance jet instabilities, leading to kiloparsec scale oscillations and filamentary structures in the jets. These filaments can appear as twisted structures that, in some projections, resemble a helical topology similar to the observed radio morphology of 3C 75. The simulations in Molnar et al. (2017) indicate that the radio morphology of the north-west side of 3C 75 may be due to a collision between two jets with distinctly different speeds and a small impact parameter. Strong instabilities generated in a fast jet colliding with a much slower jet cause the fast jet to break up into two oscillating filaments, the morphology of which resembles a double helix.
In this paper, we present a numerical study motivated by the striking morphology of the 3C 75 radio source. We investigate the evolution of two two-sided jets colliding within a stratified atmosphere in the presence of a wind representing the relative motion of the AGN as they travel through the intracluster gas. The primary aim of this work is to investigate the behaviour of merging jets, with a particular focus on how the dynamics of the ICM impacts the jet trajectories, the morphology of the jet backflow, and the coalescence of the colliding jets. The trajectories of distorted and bent jets can be used to probe the conditions in their environments. Bent and distorted jet morphologies for example can be used to reveal the substructure and dynamical motions of the ICM and shed light on the relative motions between the host galaxies and the ICM.
In our work, we build upon the numerical study of jet-jet collisions presented in Molnar et al. (2017). We extend the parameter space of the investigation presented in Molnar et al. (2017) to include a non-uniform ICM with ICM winds of various velocities and inclination angles, alongside different impact parameters of the colliding jets. We also perform the simulations in a larger computational volume, enabling the interactions of both ends of the source to be modelled. Finally, we relate our results to 3C 75 to provide insight into the complex nature of the source.

N U M E R I C A L S E T U P A N D S I M U L AT I O N PA R A M E T E R S
To investigate jet-jet collisions and the large-scale bending observed in 3C 75, we conduct simulations of two interacting two-sided jets propagating through a stratified ambient cluster medium in the presence of an external wind. The external wind represents the relative motion between the source and surrounding intracluster gas, which causes the swept-back morphology of the jets on kiloparsec scales (see Fig. 1).
We use the FLASH (magneto-)hydrodynamics code (Fryxell et al. 2000) in pure hydrodynamic mode to conduct the simulations. The simulations neglect the effects of magnetic fields, viscosity, and radiative cooling. We do not include the latter in our simulations as the cooling time in the central regions of the cluster is significantly longer than the evolutionary period that we simulate. Furthermore, we assume that the large-scale jet dynamics are governed by hydrodynamics.
The FLASH code is a Eulerian finite-volume, Godunov based, parallel code that utilizes parallel adaptive mesh refinement (AMR). FLASH supports a variety of hydrodynamic numerical schemes and in this work, we adopt the HLLC (Harten-Lax-Van-Leer with Contact) Riemann solver (Li 2005) and piecewise parabolic data reconstruction (PPM, Colella & Woodward 1984), the latter providing third-order spacial accuracy. The simulations are conducted in three dimensions in a Cartesian coordinate system.
In the absence of magnetic fields, cooling/heating and viscosity FLASH solves the Euler equations. In a Cartesian geometry, the Euler equations can be expressed in conservative form as where p and ρ are the gas pressure and density, respectively. g is the gravitational acceleration, which can be calculated from an external gravitational potential through g = −∇ . The total energy per unit mass E is the sum of the internal energy per unit mass and kinetic energy per unit mass and is given by E = + (1/2)v 2 . The kinetic energy per unit mass is given by (1/2)v 2 = 1/2(u 2 + v 2 + w 2 ), where u, v, and w are the fluid velocities in the x, y, and z directions, respectively. Equations (1), (2), and (3) represent the conservation of mass, momentum, and energy, respectively. The set of equations is closed by using the equation of state for an idea gas where γ is the ratio of specific heats. In this work, we adopt γ = 5/3, in line with Molnar et al. (2017), for both the jet material and the ambient medium. The simulation domain corresponds to a cube with dimensions of 358.36 kpc on each a side. The grid uses eight levels of adaptive refinement (AMR), and the most highly refined grid cells have dimensions of 0.175 kpc a side. We impose outflow boundary Table 1. Parameters for the different ICM wind models used in the simulations. The simulations are divided into four groups as given by the Type column, where each type represents the parameter being changed between the simulations in each set. The wind velocity, v w , is the relative velocity between source and ICM. For simulations in the wind angle group, the inclination angle of the wind is changed while the wind velocity and the impact parameter of the jet collision remain fixed. θ w is the wind inclination angle with respect to the positive x-axis (see Fig. 2). P is the impact parameter of the jet collision expressed in units of the jet diameter. The single jet simulations have been included as a control so that we can isolate the effects of the colliding jets from those of the wind.

Type
Run To explore the effects of the external wind on the interacting jets, we perform a number of simulations in which we vary the direction and velocity of the external wind, as well as simulations in which we vary the impact parameter between the colliding jets. The simulations are summarized in Table 1. The names of the simulations are formatted as follows: RT<wind inclination angle θ w > V<wind velocity v w > P<impact parameter> where the wind velocity in the simulation name is given in units of 10 km s −1 . The impact parameter in the simulation name is given in units of the nozzle diameter and is measured from the nozzle centres. The simulations comprise three different groups as specified in the first column in Table 1. In the 'wind angle' group, we change the inclination angle of the wind θ w , while fixing the wind velocity to v w = 1200 km s −1 . The inclination angles of the wind for each simulation are given in the fourth column of Table 1. We take as our zero angle for the wind inclination the vector with origin at the grid centre and aligned in the positive x direction as shown by the angle θ w in the schematic diagram in Fig. 2. In the 'wind velocity' group of simulations, we fix the wind inclination angle to 130 • and vary the wind velocity. In the 'impact parameter' group, we fix the wind inclination angle and velocity to 130 • and 1200 km s −1 , respectively, and vary the impact parameter P between the colliding jets. These three sets of simulations allow us to explore how the evolution of interacting/colliding jets depends The centre of the cluster is placed at the grid origin shown by the black 'X'. The point midway between the two jet nozzles is offset from the cluster centre by a radial distance of 48.9 kpc. The velocity vectors of the jets on the right of the simulation domain are aligned such that there is a 36 • angle, α, between them. The ICM wind has a velocity v w and an inclination angle θ w as shown in grey, the specific values for each run are given in Table 1. For the simulations in which only a single two-sided jet is simulated (simulations that comprise the 'single jets' group in Table 1), only one jet nozzle is present in the simulation domain. on the magnitude and direction of the ICM wind and the distance between the jets.

Ambient medium parameters
The cluster atmosphere into which the jets propagate is described using a density distribution given by an isothermal beta model, where r = x 2 + y 2 + z 2 , r c is the core radius, ρ 0 is the central density, and the radial profile is characterized by the β parameter. The atmosphere is centred at the grid origin (0,0,0), see Fig. 2. The β, core radius, and core density parameters are presented in Table 2 and are taken from the single spherically symmetric best-fitting model to Abell 400 given in table 5 of Beers et al. (1992). The model presented in Beers et al. (1992) uses a cluster temperature of 2.5 +0.8 −0.6 keV determined from Einstein observations of Abell 400. For our simulations, we use a cluster temperature of 2.28 keV, which corresponds to the best-fitting cluster temperature obtained from Chandra observations of Abell 400 given by Hudson et al. (2006).
The pressure of the ambient medium is calculated using where μ = 0.61. The atmosphere is supported by a gravitational potential from an assumed dark matter halo. The gravitational acceleration in the radial direction due to the potential is g r = −∇ , the magnitude of which is given by We conduct a simulation of the ICM in the absence of both the jets and external wind in order to test the stability of the atmosphere. We find that the imposed gravitational potential is sufficient to keep the atmosphere stable beyond the maximum evolution period of 42 Myr that we simulate.
The relative velocity between the nuclei and the cluster is given in Hudson et al. (2006) as 1200 km s −1 and is derived using the bestfitting overall cluster temperature (2.28 keV) and a cluster sound speed of 850 km s −1 . As we adopt a value of μ = 0.61 in this work, the sound speed in the ambient medium simulated is smaller than this (773 km s −1 ; see Table 2).
A schematic diagram of the grid setup used in these simulations is shown in Fig. 2. The dimensions and resolution of the simulation domain outlined at the beginning of Section 2 are selected from the requirements that the grid is large enough that the jets do not propagate past the boundaries throughout the period of evolution, and that each jet nozzle (blue rectangles in Fig. 2) fits into one grid block.
In this work, we are interested in the extragalactic regime of jet evolution, with a specific focus on how binary AGN jets interact with the ICM and how they interact with each other. As a result, we do not model the gas distributions of the host galaxies, or the jet propagation through the galactic gas in our simulations. Simulating the jet launch and propagation through both the galactic and intracluster gas would require such a large dynamic range, in both space and time, that the simulations would not be feasible with today's computers.

Jet parameters
The full jet parameters are given in Table 3. Jet material is injected into the simulation domain through cylindrical inflow regions or 'nozzles'. The top t and bottom b designations of the two-sided jets in Table 3 correspond to the relative position of their nozzles with respect to the y-axis, as shown in the schematic diagram in Fig. 2. For each two-sided jet, we adopt a nozzle radius of 1.05 kpc. The grid resolution is such that the nozzle radius is composed of six maximally refined grid cells, where each maximally refined grid cell has sides of 0.175 kpc. The nozzle height (this height refers to the entire cylindrical inflow region containing both sides of opposing flow velocity for each of the two-sided jets) is 0.7 kpc and is composed of four maximally refined cells. Table 3. Jet parameters in physical units for the two bipolar jets being simulated. Jet t is the weaker top jet, and Jet b is the stronger bottom jet (top and bottom are referenced from their y position on the grid). The nuclei of 3C 75 are offset from the X-ray emission-weighted centre of the cluster by 1.68 arcmin (≈48.9 kpc; Hudson et al. 2006) and the nuclei are separated by 8 kpc (deprojected separation; Hudson et al. 2006). The beta model cluster atmosphere is centred at the grid origin.

Parameter
Jet Temperature (keV) 113.24 113.24 Temperature (K) 1.31 × 10 9 1.31 × 10 9 Number density n 0 (cm −3 ) 3 . 0 × 10 −5 3.0 × 10 −5 Density (g cm −3 ) 3 . 0 6 × 10 −29 3.06 × 10 −29 Pressure (dyne cm −2 ) 5 . 4 4 × 10 −12 5.44 × 10 −12 Sound speed c j (cm s −1 ) 5 . 4 4 × 10 8 5.44 × 10 8 Velocity (cm s −1 ) 3 . 8 9 × 10 9 (0.13c) 7.00× 10 9 (0.23c) Mach number a 7.15 12.87 γ 5/3 5/3 Bipolar kinetic power b (erg s −1 ) 5 . 9 5× 10  The two two-sided jets in the simulations have a total combined kinetic power of P j = 4.06 × 10 44 erg s −1 . We derive P j using the flux densities of 3C 75 observed at 1.4 GHz presented in Ledlow & Owen (1995) and the scaling relation between the radio power at 1.4 GHz to jet kinetic power given in Cavagnolo et al. (2010). We emphasize that the jet powers derived in this way are taken purely as a starting point for simulating 3C 75, and that other scaling relations between the radio power observed at a given frequency and the kinetic power of the jet can similarly be used to provide an estimate of the kinetic jet power (e.g. Willott et al. 1999;Bîrzan et al. 2004;O'Sullivan et al. 2011;Kokotanekov et al. 2017). It is important to emphasize that the mechanical powers obtained from the observed radio powers using these scaling relations are estimates of the true mechanical power of the source. Since 3C 75 consists of two bipolar jets that interact significantly, and cocoons that appear to merge, it is not clear how (or if) such interactions affect the relationship between the radio power and the cavity power of the source and thus the derived kinetic power. The true mechanical power of 3C 75 may differ significantly from the values presented, and realistically as a result of the inherent scatter in the scaling relations and the significant interaction between the jets in the source, the mechanical jet power P j presented here is estimated at best to within an order of magnitude.
In our previous paper (Molnar et al. 2017), we found that after two colliding jets 3 merge, one jet can take over the trajectory of the other only if the former is significantly more powerful. In this case, the merged jet will have a large bending angle relative to the weaker jet, as seen in the radio image of 3C 75. The ratio between the velocities of the two bipolar jets, and thus the ratio between the kinetic powers of the bipolar jets, is taken from the simulation presented in Molnar et al. (2017) that most closely resembles the morphology inferred from the VLA radio maps of 3C 75. Correspondingly, for the simulations presented in this work, the jet density is held constant, while the jet velocity is varied between the bipolar jets (see Table 3) in order to unevenly distribute the total kinetic power of the source between the two jet sources. The bottom, powerful, bipolar jet is initialized with a velocity 1.8 times that of the weaker (top) bipolar jet, corresponding to conditions closely matching those presented in Molnar et al. (2017).
In line with the scenario observed for the nuclei of 3C 75, the jet nozzles are separated by a (deprojected) distance of 8 kpc, and we offset the jet nozzles from the centre of the cluster atmosphere. The Chandra observations presented in Hudson et al. (2006) show a projected offset between the nuclei of 3C 75 and the diffuse X-ray emission-weighted centre of the cluster of 1.68 arcmin. Since the projection angles of the nuclei are unknown, we use the projected offset to displace the jets from the centre of the isothermal-beta atmosphere in our simulations. The nozzle positions are given in the last row of Table 3. The angle α between the colliding jets (see schematic in Fig. 2) is 36 • as in Molnar et al. (2017). This is consistent with the angle between the jets in 3C 75 inferred from the VLA radio maps (e.g. see Fig. 1).
The detailed simulation setup (jet orientation, ambient medium parameters) is relevant to 3C 75. However, as we vary the wind direction across a large region of parameter space and vary both the wind velocity and the impact parameter of the jet collision, the scenarios simulated may be applicable to other jet-jet collisions/interactions resulting from cluster mergers.

Jet fractions and cocoon definition
In these simulations, we assign different passive fluid tracers to each bipolar jet. This allows us to trace how the injected jet plasma mixes during the simulation. The volume fractions associated with each jet tracer are given by f j, t and f j, b , where the subscripts t and b denote the top and bottom jet (with respect to the y-axis) in the simulation domain. Since two jet tracers are used in these simulations, the numerical definition of the cocoon and the combined bulk jet flow of the two bipolar jets in this work are based on the volume fractions of the ICM material f ICM . The cocoon is therefore identified as the volume for which f ICM < 0.999, as shown by cyan curves in Fig. 3. The regions that we identify as the bulk jet flow of both bipolar jets is defined as the volume inside of which f ICM < 0.2 (black curves in Fig. 3), where a non-zero value is chosen to take account of any mixing between the jet and the ICM. We note that these definitions are arbitrary and will be impacted by the level of numerical diffusion present in the simulations. However, we find that they identify their respective regions well across the entire evolution, compared to other values of f ICM that we have considered.

S I N G L E J E T E VO L U T I O N
First, we outline the jet evolution in simulations RT130 V120 top and RT130 V120 bot that comprise the 'single jet' simulations shown in Table 1. These simulations evolve the powerful bipolar jet (run RT130 V120 bot, with a total jet power P j = 3.47 × 10 44 erg s −1 ) and the comparatively weaker bipolar jet (run RT130 V120 top, with P j = 5.95 × 10 43 erg s −1 ) separately within the isothermal beta atmosphere. Each bipolar jet propagates through the ICM in the presence of a cross wind with velocity of 1200 km s −1 and an inclination angle of 130 • .
In RT130 V120 top and RT130 V120 bot, the supersonic advance of the jets into the ambient medium drives a pair of shocks: a forward-moving bow shock that precedes the jets and a reverse shock that acts to slow the jet advance into the ambient medium. Jet material that crosses the reverse shock flows back along the The cocoon region is given by the volume inside of which the ICM tracer f ICM satisfies f ICM < 0.999 as shown by the cyan contour. The black contour shows f ICM < 0.2, which is used to identify the bulk jet flow. The location of each jet nozzle is shown by a black X. The colour scale shows the log density between 1.52 × 10 −30 g cm −3 and 2.28 × 10 −27 g cm −3 . Regions of the jets furthest away from the nozzles are highly disrupted and can frequently move out of the central plane, truncating the black contours -this is particularly the case for run RT130 V120 P0 as shown in panel b. jet inflating a cocoon of shocked jet material. The bow shock compresses the ICM that crosses it, generating a shocked shell of ambient medium around the cocoon backflow and bipolar jet. A contact discontinuity separates the shocked shell from the cocoon of shocked jet material. Due to the velocity shear between the shocked ICM shell and the cocoon, the contact discontinuity is susceptible to Kelvin-Helmholtz instabilities, which cause shocked ambient medium to mix with the shocked jet material in the cocoon.
During the early evolution of the jets (t ≤5 Myr), the supersonic expansion of the cocoon, driven by the high internal cocoon pressure, dominates the ICM wind and the cocoon expands symmetrically around the jets. In this early evolutionary phase, the jets remain relatively straight as they propagate and terminate abruptly at the jet shock. This initial phase of the jet evolution is consistent with that of a jet propagating through a uniform or stratified ambient medium in the absence of an ICM wind, as found in the hydrodynamic and magneto-hydrodynamic simulations conducted by Reynolds, Heinz & Begelman (2002), Zanni et al. (2003) (2019), and references therein.
The ICM wind in runs RT130 V120 top and RT130 V120 bot has a large perpendicular velocity component with respect to the initial (unperturbed) direction of jet propagation. The jets are therefore subjected to a transverse ram pressure that remains unbalanced by the ram pressure of the jet over the course of the evolution. The presence of the wind will therefore (1) significantly modify the morphology of the cocoon and (2) induce large-scale curvature in the jets.
As the evolution progresses, the internal pressure of the cocoon approaches pressure equilibrium with the surrounding ICM. The initially supersonic expansion of the cocoon slows and becomes insufficient to counteract the ram pressure of the wind. The cocoon on the side of the jet exposed to the wind is crushed. Jet backflow displaced from the wind-facing side of the cocoon then accumulates on one side of the source. Cocoon material is increasingly deflected by the wind as the evolution time increases and is preferentially dragged in the direction of the wind, creating turbulent plumes or tails as shown in Fig. 4.
Regions of the jet that are pushed against the inner edge of the cocoon drive a shock into the surrounding shell of shocked ICM. As a result of this shock, the jet flow on the wind-facing side of the jet is confined by a high-pressure gradient that redirects the jet trajectory. Confinement by the pressure gradient leads to the largescale bending of the jets on either side of the jet nozzle. By t ≈ 10 Myr, the jets in both RT130 V120 top and RT130 V120 bot have started to bend along their axes. The pressure gradient also inhibits the lateral expansion of the bulk jet flow, resulting in enhanced collimation of the jet on the side of the cocoon exposed to the wind.
Multiple oblique shocks are driven into the jet following the onset of jet bending. The shocks destabilize the jet, promote the growth of instabilities within it, and act to destroy its collimation. Beyond the point of maximum curvature along the jet axis, the jet flow becomes more turbulent as Kelvin-Helmholtz instabilities are further excited. As a result, the jet flow is increasingly disrupted along its length and broadens, spreading its momentum flux over a larger area and thus slowing the jet's advance into the ambient medium. Consequently, the jet shock weakens with increasing evolution time until the jets no longer terminate abruptly but rather become increasingly turbulent and gradually transition into the diffuse (with respect to the value  Table 3 for parameters). In both runs, the bipolar jet is simulated the presence of a wind with a velocity of 1200 km s −1 and an inclination angle of 130 • . of the jet tracer fraction, see dark orange material in Fig. 4) tails or plumes of jet material that widen downstream of the jet.
The sound speed (and temperature) of the jet material is highest in the plume regions at the ends of the jets. The plumes therefore respond the fastest to a disturbance by the wind. As the evolution time increases, the leading end of the jets (and much of the plume material) is preferentially swept in the direction of the wind. This accelerates material within these flow regions, since their propagation direction becomes roughly aligned with the direction of the wind. Consequently, jet material accumulates at the end of the plumes as shown by the purple and blue regions in Fig. 4.
The disruption of the jet downstream of the location of maximum curvature along the jet axis can cause the jet to break apart into filamentary structures, which are then lifted from the bulk jet flow by the ICM wind. An example of this behaviour can be seen in Fig. 4 (top left-hand panel) for the jet propagating to the left in run RT130 V120 top. Despite the onset of turbulence and disruption of the bulk jet flow downstream of the nozzles, the jets in simulations RT130 V120 top and RT130 V120 bot typically remain coherent enough to survive beyond propagation distances of ∼50-120 kpc following their point of maximum deflection.
Due to the continual disruption of the jets downstream of the nozzles and the subsequent absence of a strong jet termination shock, the jets do not display hotspots -regions of localized high temperature at the head of the jet, bracketed by the bow and jet shocks that are common features of powerful radio sources. Hotspots correspond to compact sites of strong radio emission. The jets in WAT sources, such as 3C 75, are thought to terminate in a variety of ways; an initially well collimated WAT jet may display hotspots at the outer edges of the plume at locations where the jet happens to impinge on the plumes edge, or the jets may not display hotspots if (1) they smoothly transition into the plume or (2) propagate without disruption for some distance into the plume before 'disappearing' inside it (Hardcastle & Sakelliou 2004). The absence of hotspots in WATs is not a feature of all observed WATs: some display compact hotspot features (e.g. see Hardcastle 1998).
Briefly comparing the evolution in runs RT130 V120 top and RT130 V120 bot with previous works in which bipolar jets are simulated in the absence of an external wind (e.g. Reynolds et al. 2002;O'Neill et al. 2005;Gaibler et al. 2009;English et al. 2016), we find that: (1) the wind perturbs the jets to such an extent that the jet shock weakens and disappears following the early evolutionary period.
(2) For jets of comparable kinetic power, the leading ends of the jets are more highly disturbed when a wind is simulated. (3) The wind enhances the onset of instabilities within both the bulk jet flow and at the shear layer between the cocoon and shocked ICM. This will likely promote mixing between the jet and the ICM material in simulations in which a wind is simulated. (4) A sequence of shocks is driven into the jets at the location of maximum curvature. (5) The wind deflects the slightly lower velocity material at the outer regions of the bulk flow. This can distort the profile of the jet beam and lead to filamentation of the flow.
We note that the jets simulated are purely hydrodynamic in nature, which could affect the way in which they respond to the wind compared to those that are magnetized. For example, a magnetic field of suitable strength and coherence could act to stabilize the jet against instabilities (e.g. Appl & Camenzind 1992), reducing its decollimation following the deflection by the wind, which may in turn reduce the degree of jet curvature present in the simulations. However, a strong toroidal component of the magnetic field on the other hand may destabilize the jets by generating magnetohydrodynamic instabilities within it (Bodo et al. 1994(Bodo et al. , 1995(Bodo et al. , 1998(Bodo et al. , 2019. These instabilities can destroy the jets over time.

C O L L I D I N G J E T E VO L U T I O N : N O I C M W I N D
We restrict the following discussion to simulation RT000 V000 P0 in which colliding jets propagate in the absence of any ICM wind. In this run, the impact parameter P = 0 (no offset between the jets along the z-axis) and the jets collide head-on. The advance of dual two-sided jets propagating in close proximity to one another, with converging trajectories, significantly changes the evolution and structure of the jets, shocked shell, cocoon, contact discontinuity, and bow shocks compared to the runs in which only one two-sided jet is simulated.
At early evolution times, the bow shocks driven into the ambient medium by the jets propagating to the right of the simulation domain merge, due to the converging jet trajectories and impact parameter of P = 0. The jets travelling towards the right of the domain then propagate within a single bow shock and a common cocoon. The shocked backflow from both jets mixes within the common cocoon as shown in panel (a) of Figs 3 and 5.
By t 1 Myr, the bulk flow of the jets propagating towards the right of the domain has collided. Immediately following the collision, the jets expand where they meet, as the collision disrupts the collimation of the jet flow (see Section 4.1). The jets (bulk flow) are susceptible to Kelvin-Helmholtz instabilities, which can result in the distortion of the jet surface (surface modes), or the entire beam (reflected body modes; Birkinshaw 1991). The surface modes of the Kelvin-Helmholtz instability decrease in amplitude sharply as a function of distance from the jet surface. Reflected body modes dominate the entire gas in the jet and dominate if the jet Mach number M j 2 √ 2 (e.g. Ferrari 1998). The jets simulated have internal Mach numbers M j of 7.2 (top, weak jet) and 12.9 (bottom, powerful jet), and so we expect the reflected body modes to dominate. The instabilities can result in large-scale disruption of the jet by inducing shocks within the flow and exciting asymmetric and helical perturbations (e.g. Bodo et al. 1994, and references therein), which can eventually lead to its complete destruction. Again, the inclusion of strong magnetic fields is expected to suppress the growth of the Kelvin-Helmholtz instabilities and could therefore have a significant impact on the stability of the jets following their collision.   The cocoons and bow shocks associated with the jets propagating towards the left of the simulation domain do not fully merge, and two cocoons can clearly be identified downstream of the jet nozzles (e.g. see Fig. 3, top left panel of Fig. 6 and panels (e)-(h) of Fig. 7). A section of shocked ambient medium lies between the two cocoons and expands with increasing evolution time.

Collision region
In the following discussion, we define the region comprising the combined bulk jet flow of both two-sided jets as the volume in which the ICM tracer f ICM < 0.2. (The value is non-zero to allow for mixing between the injected jet material and ICM material entrained into the cocoon.) The collision region is then taken to be the point at which the bulk jet flow of the top and bottom jets propagating towards the right of the simulation domain collides. We identify the point of collision by computing the location of maximum pressure that is both within the cocoon and within a radius of 10 kpc from the mid-point between the jet nozzles.
For all simulations with zero-impact parameter (P = 0, both with and without the ICM wind, and for all wind parameters), the collision of the jets propagating to the right of the simulation domain is violent and the bulk jet flows merge together following the point of collision. The collision produces a standing shock at the interaction point between the bulk flow of the jets, which persists for the duration of the evolution. The standing shock is characterized by a local increase in density and also pressure (exceeding five times the ambient pressure and typically largest for jets propagating through an approximate cross wind) as shown by the white regions in the jets at the collision point in Fig. 8. This local increase in pressure and density would correspond to an enhancement of the magnetic field and may act as a site for possible particle acceleration and enhanced jet brightness. The standing shock, and subsequent shocks along the jet axis, can be seen in Fig. 9, which displays the ratio of the magnitude of the pressure gradient to the pressure. This ratio can be used as a proxy for shock strength.
The jet collision is accompanied by a loss of stability, collimation, and a quick transition into a highly turbulent jet flow. The generated turbulence is strong, resulting in the heating of the post-collision jet flow downstream of the collision point, alongside enhanced instabilities along the jet axis. Consequently, the resulting 'merged' jet following the collision is significantly wider than the two incident jets. The instabilities excited in the merged jet become nonlinear and increase in amplitude with evolution time and distance along the jet axis. As a result, the bulk jet flow following the collision displays kiloparsec-scale oscillations, kinks, wiggles, and filamentary structures along the jet axis. These features are transient and are most prominent in regions of the jet flow that are furthest away from the nozzles.
We find that the oscillations of the bulk jet flow following the collision (merged jet) are three-dimensional and are not the result of an in-plane instability, as the plane of oscillation is found to change during the evolution. Additionally, any waves induced in the merged jet following the collision are not rotating, and so the oscillations are not the result of a chiral instability. Thus, the bulk jet flow following the collision oscillates randomly, 'whipping' or 'flapping' around in the cocoon. Due to the oscillating jet trajectory, Figure 8. Two-dimensional pressure slice through the jet centres at t = 38 Myr. The panels correspond to (a) run RT000 V000 P0 (no wind), (b) run RT130 V120 P0 (θ w = 130 • , v w = 1200 km s −1 ), (c) run RT225 V120 P0 (θ w = 225 • , v w = 1200 km s −1 ), and (d) run RT315 V120 P0 (θ w = 315 • , v w = 1200 km s −1 ). Colour scale shows log pressure between values of 1.1 × 10 −13 dyn cm −2 and 2.2 × 10 −11 dyn cm −2 in each image. White corresponds to regions with high pressure and dark orange/black corresponds to regions with low pressure. Note the localized regions of high pressure within the shocked shell of ICM. jets can bounce off, or impinge on, the inner edge of the cocoon for simulations in which the wind is roughly perpendicular to the initial unperturbed trajectory of the jets. Regions of the jet that make contact with the inner cocoon edge display a local temporary increase in jet pressure. This increase in jet pressure would lead to an increase of the radio emissivity of the jet plasma in these regions and thus transient compact regions of enhanced radio emission, or 'knots', may be observed at the locations where the jet comes into contact with the inner cocoon edge. This behaviour is in line with the idea presented by Hardcastle (1999) in which WAT jets flapping within the plumes show only compact termination features when they happen to impinge on the edge of the cocoon/plume (see also Hardcastle & Sakelliou 2004).
Following the collision and the resulting standing shock, subsequent shocks driven along the axis of the merged jet are weaker, reducing in intensity further away from the collision point as shown in Fig. 9. The resulting trajectory of the merged jet following the collision is similar to the original trajectory of the powerful jet prior to the collision. This result is consistent as with Molnar et al. (2017). Again, none of the jets simulated (on both sides of the simulation domain) terminate in the strong jet shock typical of powerful FRII radio jets.

Changing wind angle
The simulations in which the wind angle is varied fall under the 'wind angle' group shown in Table 1 (excluding run RT000 V000 P0). Simulations in which the wind angle is roughly parallel or antiarallel to the jet trajectory at initialization (head/tail winds, e.g. RT045 V120 P0 and RT225 V120 P0) evolve similarly to the run conducted without the wind, though displaying different localized regions of high pressure (and density) in the surrounding shocked shell of ICM (e.g. see Fig. 8). These regions will correspond to locations of enhanced X-ray surface brightness within the shocked ICM shell and should lead to interesting X-ray structures surrounding the cavity inflated by the jets.
Winds that are roughly perpendicular to the instantaneous jet direction (cross winds, e.g. runs RT100 V120 P0, RT130 V120 P0, RT150 V120 P0, RT180 V120 P0, and RT270 V120 P0) have a considerable impact on the degree of jet curvature and stability, the degree of mixing between the jet and ambient fluid, and the large-scale morphology of the jets, plumes, and shocked ICM shell.
In Fig. 7, we show the morphology of the jet, cocoon region, and shocked shell of ICM at various evolution times for runs RT130 V120 P0 and RT000 V000 P0. At early evolution times during the initial supersonic expansion of the overpressured cocoon [t = 2.5 Myr, panels (a) and (e) of Fig. 7], differences between the jet morphology in runs RT130 V120 P0 and RT000 V000 P0 are minimal. As the cocoon expansion slows, the component of the wind perpendicular to the jets in run RT130 V120 P0 crushes the cocoon on the windward side of the source as previously outlined for the single jet runs. By t ∼ 10 Myr [panel (b) of Fig. 7], the cocoons on both sides of the simulation domain have merged and the jets share one common cocoon for the remainder of the evolution. Cocoon material from the jet directly exposed to the wind is increasingly deflected into the path of the weak jet by the wind. This leads to significant mixing of the cocoon backflows on both sides of the source, as shown by the passive jet tracers in Fig. 5 and density slices in Fig. 7 (top row), affecting the evolution and appearance of the weak jet.

Changing wind velocity
In the 'wind-velocity' group in Table 1, we vary the relative velocity between the jets and the ICM from 480 km s −1 to 3000 km s −1 . The wind angle fixed to 130 • for all simulations in this group. Increasing the wind velocity amplifies the degree of jet disruption and curvature (see the first column of Fig. 6). The jets are deflected, and transition into plumes at shorter propagation distances from the nozzles as the ram pressure of the wind is increased. Plume material is more effectively swept in the direction of the wind, and stronger shocks driven into the bulk jet flow as the wind velocity is increased. The latter can be seen in Fig. 9 where run RT130 V300 P0, in which the largest wind velocity is simulated, displays the strongest shocks within the jet (see red and orange regions in Fig. 9). The jets in run RT130 V300 P0 are deflected by almost 90 • from their initial propagation direction and the source morphology begins to resemble that of a Narrow Angle Tail (NAT) jet. NATs are usually associated with galaxies moving at high velocities in the gravitational potential of the cluster (Miley 1980). Jones et al. (2017) conduct 3D magneto-hydrodynamic simulations of NAT radio jets. The simulations presented in Jones et al. (2017) show that NAT jets remain fairly coherent once they are bent (unlike the morphology we obtain in run RT130 V300 P0) but can generate patchy magnetic filaments in their tails. Additionally, the authors note that the jet tails can also flap around -we find that this behaviour occurs in all our simulations in which a cross wind is present.
Large-scale oscillations of the bulk flow are amplified and jet filamentation is enhanced as the wind velocity is increased. For the largest wind velocities simulated (v w = 1200 km s −1 and 3000 km s −1 ), the wind can lift the filaments off the bulk jet flow, deflecting them in the direction of the wind. Any oscillations in these filaments may then cause the filaments to cross, or come spatially close to, the bulk flow of the jet. When viewed in some projections, this behaviour can lead to regions with high jet tracer fractions interspersed with regions of low jet tracer fractions. This would correspond to regions of bright radio-emitting jet plasma interspersed with regions of comparatively low surface brightness radio emission, or regions of Faraday rotation 'banding': Rotation measure gradients between the centres and edges of the jet profile, resulting in a 'banded' pattern of the rotation measure (e.g. Laing & Bridle 1987).

Changing impact parameter
We now consider the simulations in which we vary the impact parameter P between the jets ('impact parameter' group in Table 1). The wind velocity and inclination angle are fixed to 130 • and 1200 km s −1 , respectively, for all simulations in the 'impact parameter' group. Due to resolution constraints, the smallest non-zero impact parameter we can simulate is P = 1 (run RT130 V120 P1), which corresponds to an offset between the nozzle centres in the z direction of one jet diameter (2.1 kpc). Fig. 10 shows the projections of the ICM tracer f ICM for simulations RT130 V120 P0, RT130 V120 P1, and RT130 V120 P3 (P = 0, 1, and 3, respectively). The bulk jet flows on the right of the simulation domain in simulations RT130 V120 P1 and RT130 V120 P3 do not make contact as the wind is planar, and the jets do not deviate significantly from their initial trajectory in the z direction over the course of the evolution. Due to the absence of a jet-jet collision in runs RT130 V120 P1 and RT130 V120 P3, the entire source resembles an 'X'-shaped structure as the weak and powerful jets do not bend to the same degree. Consequently, the jets advance into the ambient medium with two interacting bow shocks on either side of the simulation domain.
Although there is no direct contact between the bulk jet flows in runs RT130 V120 P1 and RT130 V120 P3, their respective cocoons still make contact and the jets interact through turbulence: Turbulent feedback from the interacting cocoons increases the level of jet disruption compared to the simulations in which only a single two-sided jet is simulated. The effect of turbulent cocoon feedback on jet structure is most evident in run RT130 V120 P3 (see top right-hand panel of Fig. 10) where the weak jet breaks apart into kiloparsec-scale filaments following the point of maximum jet deflection.

Implications for 3C 75
The single jet runs RT130 V120 top, RT130 V120 bot and the impact parameter runs RT130 V120 P1, RT130 V120 P3 alongside run RT130 V120 P0 are all simulated with the same wind parameters (v w = 1200 km s −1 , θ w = 130 • ). A comparison between these runs enables us to determine the impact of the close proximity of another jet on the resulting large-scale jet morphology. We conclude that: (1) The jets still interact via turbulence even if the bulk flows do not make contact, leading to broader jets that become more unstable along the jet axis. (2) The collision amplifies jet disruption and kiloparsec-scale oscillations in the bulk flow. (3) Without some contact between the bulk jet flows, the location of maximum bending for the weak and strong jet propagating to the right of the simulation domain does not occur at the same location in the atmosphere: With the non-zero impact parameters simulated, the bulk flows of the jets will appear to cross each other in projection, leading to an 'X'-shape structure that is not observed in 3C 75. The simulations therefore imply that some contact is required between the bulk jet flows in order to produce a morphology consistent with 3C 75.
Throughout this work, we assume that both of the two-sided jets have the same density, alongside an uneven distribution of the total kinetic power of the source between them. These assumptions have considerable implications on the resulting morphology of the jets. The top jet, for example may be less dense and/or significantly weaker than the bottom jet and thus more easily deflected by the wind -this may result in the jets on the right of the simulation domain bending in tandem as observed in 3C 75, even in the absence of a collision.

Hybrid morphology?
In some projections (e.g. see the y-and z-axis projections in Columns 2 and 3 of Fig. 6), the asymmetric inflation of the cocoon and the deflected jets can appear to be a single two-sided jet with vastly different morphology on either side of the jet origin. A large disparity in the radio-jet morphology on either side of the nucleus is characteristic of HYbrid MOrphology Radio Sources (HYMORS). HYMORS exhibit different Fanaroff and Riley (FR) morphology (Fanaroff & Riley 1974) on either side of the nucleus. In HYMORS, the jet on one side of the nucleus exhibits a typical FR-I type morphology; poorly collimated jets displaying diffuse plumes of radio emission that are brightest at distances within the inner half of the radio source, while the other jet exhibits an FR-II type morphology, displaying highly collimated jets terminating in hotspots and complex lobes. HYMORS are extremely rare; of the 1700 sources in the VLA FIRST survey (Faint Images of the Radio Sky at Twentycentimetres; White et al. 1997), Gawroński et al. (2006) find that HYMORS constitute ∼ 1 per cent or less of the sources. Gopal-Krishna & Wiita (2000) suggest that HYMORS are the result of kiloparsec scale interactions between the jet and an ambient medium with dissimilar properties on either side of the jet nucleus, rather than any intrinsic asymmetric differences in the central engine. This hypothesis is supported by Chandra observations of HYMORS (e.g. Miller et al. 2006;Miller & Brandt 2009) and highresolution VLBA observations, which appear to suggest that the jets  Fig. 4. The columns correspond to simulations conducted with a wind velocity and direction of 1200 km s −1 and 130 • , respectively, and an impact parameter between the two colliding jets of P = 0 (left, run RT130 V120 P0), P = 1 (middle, run RT130 V120 P1), and P = 3 (right, run RT130 V120 P3) where P is given in units of the jet nozzle diameter. P = 0 signifies a head-on collision between the jets propagating towards the right of the simulation domain, when P = 1 the jets just miss each other, and when P = 3 the jets are significantly offset from each other. The simulations show that in the absence of a direct collision between the converging jets, the jets do not bend in tandem and simply appear to cross. The bending of the upper left jet towards the upper right-hand corner of the simulation domain is amplified with a non-zero impact parameter (see top rows) as it becomes increasingly less protected from the wind by the bottom jet.
in HYMORS on the 1-10 kpc scale resemble FR-IIs (Cegłowski, Gawroński & Kunert-Bajraszewska 2013). However, in the absence of any clear large-scale anisotropies in the ICM on either side of the source, we find that jet-jet interactions are another way of creating asymmetrical jet structures when viewed in some projections. Fig. 11 shows the velocity of the bulk flow of the merged jet following the collision (the velocity is given in units of the sound speed c s in the undisturbed ambient medium, c s = 7.72 × 10 7 cm s −1 ) at t = 30 Myr. The velocity is shown as a function of position along the jet. The position along the jet is measured from the location of the collision point in each simulation. In order to compare the deceleration of the merged flow to that of the single two-sided jet runs, we also compute the velocity along the jet axis for runs RT130 V120 top and RT130 V120 bot. For RT130 V120 top and RT130 V120 bot, we measure the position along the jet from the average x-coordinate of the jet collision in the remainder of the simulations plotted. Each of the curves displayed in Fig. 11 ends at the position where the numerical definition of the bulk jet flow (the volume in which f ICM < 0.2) ends, and thus the curves terminate at different lengths for each simulation.

P O S T-C O L L I S I O N J E T D E C E L E R AT I O N
The top and bottom bipolar jets are initialized with velocity v j /c s = 50.4 and 90.7, respectively, in each simulation. For run RT130 V120 bot in which only the powerful two-sided jet is simulated, regions of the jet that still fall under the numerical definition of the bulk flow at t = 30 Myr are hardly decelerated from the jet injection velocity of v j /c s = 90.7. Run RT130 V120 top (solid yellow curve, Fig. 11), in which only the weaker top jet is simulated, displays more variation in the velocity across the jet axis. The dip in velocity at ∼10 kpc for run RT130 V120 top corresponds to the location of maximum curvature along the jet axis. Following this location, the bulk flow starts to accelerate, as it is increasingly re-orientated in the direction of the wind downstream of the nozzles. This flow acceleration is similarly the case for run RT100 V120 P0 (solid blue curve).
For simulations in which the jets collide (solid black, red, blue, and grey curves, dashed green curve in Fig. 11), the resulting merged jet is typically continually decelerated along its length. Simulation RT225 V120 P0 (solid grey curve), in which the merged jet propagates against a head-wind, consistently displays the lowest velocity compared to the other simulations in which the jets collide. Figure 11. Velocity of the bulk flow (numerically defined as f ICM < 0.2, see Section 2.3) along the merged jet after the point of collision at t = 30 Myr for simulations RT000 V000 P0, RT045 V120 P0, RT100 V120 P0, RT130 V120 P0, and RT225 V120 P0 (see Section 6 for more details). In simulations RT130 V120 top and RT130 V120 bot, only one bipolar jet is simulated, and thus no collision takes place. For RT130 V120 top and RT130 V120 bot, we show the velocity along the single jet from a position corresponding to the average location of the collision in the simulations in which the jets do collide. The x-axis corresponds to the distance along the jet from the collision point onwards. The y-axis shows the jet velocity in units of the sound speed c s = 7.72 × 10 7 km s −1 of the undisturbed ambient medium at initialization. The maximum distance on the x-axis for each curve corresponds to the position at which the jet no longer satisfies the numerical definition of the bulk flow, due to mixing with the entrained ICM. The curves in this figure have been smoothed. Following the point of collision, the merged jet decelerates to velocities substantially lower than those of the incident jets (50.37 and 90.67 in units of c s for the top and bottom jets, respectively). The most severely decelerated (grey curve, run RT225 V120 P0) corresponds to the simulation in which the merged jet is propagating against a head-wind.
The bulk jet flow in run RT225 V120 P0 is decelerated to v j /c s ≈ 27 at the leading end of the merged jet. Run RT045 V120 P0 (dashed green curve), in which the post-collision jet propagates in a tail-wind, decelerates much more gradually along its length and reaches comparable velocities only when it has propagated an additional 30 kpc into the ambient medium. The shocked backflow comprising the plume material in comparison typically has a velocity in the order of v j /c s ∼ 10−20 in all the simulations presented in Fig. 11. Fig. 12 shows the ratio of the total cocoon mass M C to the total injected jet mass M inj as a function of evolution time. The ratio can be used as a measure of the entrainment of intracluster gas into the cocoon (see e.g. Reynolds et al. 2002). The total cocoon mass is computed by

E N T R A I N M E N T
where C is the cocoon volume as defined in Section 2.3. The entrainment rates begin to diverge at t 8 Myr and the simulations conducted with the ICM wind consistently display larger values of M C /M inj . Therefore, the wind increases the level of entrainment following the early evolution of the jets. The rate of Figure 12. Ratio of the cocoon mass, M C , to the total mass injected by the jets, M inj , which can be used as a measure of entrainment of the ICM. All simulations conducted with the wind show larger values of M C /M inj , signifying that the jets in these simulations entrain more ICM than the simulation conducted without the wind (run RT000 V000 P0, black curve). A non-zero impact parameter (dashed blue P = 1 and magenta P = 3 curves) results in an increase in the level of entrainment from very early times in the evolution (t ≈ 8 Myr). Simulations in which the wind has a large velocity component perpendicular to the initial propagation direction of the jets (e.g. solid red and purple curves) typically result in increased levels of entrainment compared to those in which the wind velocity is roughly head/tail-on to the direction of jet propagation. entrainment settles at a roughly constant value by t = 20 Myr for all simulations (compare the slopes of the curves in Fig. 12) regardless of impact parameter, wind direction and velocity, or whether a wind simulated.
Increased levels of entrainment in simulations with the wind are expected as (1) the jets interact more significantly with the cocoon edge (e.g. runs in which the wind has a large velocity component perpendicular to the direction of jet propagation: runs RT130 V120 P0, RT270 V120 P0, RT180 V120 P0), and (2) the wind increases the velocity shear between the cocoon material and the surrounding ICM (particularly in simulations with head/tail wind, e.g. RT045 V120 P0 and RT225 V120 P0). The larger velocity shear further excites Kelvin-Helmholtz instabilities at the contact discontinuity between the cocoon and shocked ICM shell, thus promoting mixing between the jet and ICM material.
For simulations in which the wind has a large velocity component perpendicular to the jet propagation direction, the interaction between the jet and contact discontinuity can be quite substantial. This is particularly true for high-velocity cross winds (v w ≥ 1200 km s −1 ) at late evolution times (t ≥ 35 Myr), where a prolonged interaction between the jet and contact discontinuity distorts the contact surface to such an extent that filaments of shocked ICM are effectively dragged into the cocoon. This behaviour is particularly prevalent for simulation RT130 V300 P0 in which the largest wind velocity is simulated. We find that the amount of entrained ICM depends more strongly on the wind velocity rather than inclination angle.
Simulations conducted with a non-zero impact parameter (runs RT130 V120 P1 and RT130 V120 P3, dashed blue and magenta curves in Fig. 12, respectively) display the highest levels of entrainment for the majority of the evolution. This is due to the absence of a cocoon merger on either side of the simulation domain, which enables all the jets (on both sides of the simulation domain) to Figure 13. The temporal evolution of the cocoon internal (top panels, left: early evolution times, right: full evolution), kinetic, and potential energy (bottom panels). Only a selection of simulations is shown for clarity. The internal energy dominates over the kinetic, and only a small contribution to the cocoon energy is in the form of potential energy. In simulations where two bipolar sources are evolved (all curves except the cyan and orange dash-dotted curves that correspond to the single jet simulations), the presence of the ICM wind typically increases the internal, kinetic, and potential energy in the cocoon, particularly for winds that have a large perpendicular velocity component with respect to the direction of jet propagation (e.g. run RT130 V120 P0, solid red curves). Increasing the wind velocity results in an increase of all components of the cocoon energy (compare dotted green, solid red, and dotted dark red curves). The effect of increasing the impact parameter on the cocoon energy is less significant in comparison to that of the wind velocity. make direct contact with the inner edge of their respective cocoons. The difference in the rate and level of entrainment is negligible between the non-zero impact parameters tested, at least until t = 30 Myr where they begin to diverge.

C O C O O N E N E R G E T I C S
We compute the temporal evolution of the cocoon energy in order to determine the impact of an external wind and two interacting cocoons on the energetics of the system. The cocoon internal energy is computed by the cocoon kinetic energy by while the gravitational potential energy of the cocoon is calculated by The cocoon region (C) over which the integrals are performed is defined as volume inside of which the passive ICM tracer f ICM < 0.999 (see Section 2.3 and cyan contour in Fig. 3). The definition of the cocoon therefore includes both the backflow region and the bulk jet flow (black contours in Fig. 3). The energies are referenced from their values on the grid at initialization, such that E = E(t) − E(0) for all energy components.
The evolution of the cocoon internal E int , kinetic E kin , and potential E pot energy is shown in Fig. 13 for a selection of simulations. The presence of the ICM wind generally leads to an increase in each of the cocoon energy components at a given evolution time. The bottom bipolar jet dominates the cocoon energetics, and runs in which two bipolar jets are simulated typically show similar trends to run RT130 V120 bot (cyan dash-dotted curve in Fig. 13), in which only the more powerful bottom bipolar jet is evolved unless the wind velocity is very large (v w = 3000 km s −1 , dotted dark red curve).
The internal energy increases roughly linearly with time and dominates over the kinetic and potential energy contribution to the cocoon energy for all simulations. These results are consistent with the energetics analysis of single jets propagating in the absence of an external wind presented in e.g. Reynolds et al. (2002); Vernaleo & Reynolds (2007); Gaibler et al. (2009);and Hardcastle & Krause (2013). The decay of the internal energy at late times (t ≥ 26 Myr) for simulation RT130 V300 P0 (dotted dark red curve in top left-hand panel of Fig. 13) is due to regions of the cocoon exiting the outer boundary of the simulation domain beyond these times. Regrading each energy component, all curves follow the same general trends; however, those for RT130 V300 P0 (in which the largest wind velocity is simulated) remain substantially offset. Run RT130 V300 P0 typically displays significantly higher internal, kinetic, and potential energy throughout the entire evolutionary period.
The cocoon internal energy in all simulations with two bipolar jets (both with and without a collision) evolves similarly and in line with the combined internal energy of the single jet simulations RT130 V120 top and RT130 V120 bot until t ∼ 3 Myr, despite the jets typically colliding by t 1 Myr (in the simulations in which a collision takes place). For 3 ≤ t ≤ 17 Myr, the internal energy of all simulations with dual bipolar jets grows at a faster rate than that of the combined internal energy of the single jet simulations (see black dashed-dotted curve in top right-hand panel of Fig. 13). This results in larger values of the cocoon internal energy for simulations in which two bipolar jets are simulated, compared to the total combined internal energy of the two single jet simulations. The cocoon internal energy is therefore enhanced when two jets propagate in close proximity to one another, regardless of whether their bulk flows make contact. Furthermore, all simulations in which a jet-jet collision occurs have energies larger than simulations in which the jets miss each other (P = 1 and 3) for t ≤ 17 Myr. This implies a contribution to the internal energy from both the jet collision and the interactions between the cocoons.
After t ≥ 17 Myr, the simulations that initially had the largest cocoon internal energy prior to this evolution time become those with the least for all simulations with dual bipolar jets except the outlier, run RT130 V300 P0, which continues to display the largest internal energy.
For the cocoon kinetic and potential energy, the simulations that initially display the largest energy (e.g. RT130 V300 P0, RT130 V120 P1, RT130 V120 P3, and RT130 V120 P0) continue to do so for the remainder of the evolution. The cocoon kinetic and potential energy increases with increasing wind velocity (see dotted green curve for v w = 480 km s −1 , solid red curve for v w = 1200 km s −1 , and dotted dark red curve for v w = 3000 km s −1 ). The cocoon kinetic and potential energy components are always larger for simulations with the wind compared to the simulation conducted without it and are typically comparable to (kinetic) or marginally larger still (potential) for runs with a non-zero impact parameter.

C O C O O N PA R A M E T E R S
Fig. 14 shows the cocoon temperature T c and specific entropy s c = P c /ρ 5/3 c (where P c and ρ c are the cocoon pressure and density, respectively) averaged over the cocoon volume. The quantities presented in Fig. 14 are normalized by their values at the jet injection.
We compute the volume averaged cocoon density and pressure (not shown here) and find that the evolution of these parameters does not have a strong dependence on the wind velocity, angle, and impact parameter between the jets.
The volume averaged cocoon temperature and entropy, however, show a strong dependence on the wind velocity once the jets have evolved for longer than t ∼10 Myr. Increasing the velocity of the wind (compare Fig. 14 dotted green, solid red, and dotted dark red curves for wind velocities of 480, 1200, and 3000 km s −1 , respectively) decreases the volume average cocoon temperature and entropy by as much as 64 per cent and 87 per cent, respectively, by the end of the evolution period. The fastest wind (3000 km s −1 , RT130 V300 P0 dotted dark red curve) results in cocoon material that is cooler, more dense, and of lower specific entropy. The hottest cocoon region with the highest entropy corresponds to simulation RT000 V000 P0 (solid black curve in Fig. 9) in which no wind is simulated. We therefore find that presence of wind, regardless of its inclination angle or velocity for the parameter space explored, acts to reduce the volume averaged temperature and entropy of the cocoon.
Regarding the wind inclination angle, significant changes in the volume averaged temperature and entropy can be obtained by orientating the wind such that it is almost head-on to the jet propagation direction on one side of the source. This brings the evolution of the volume averaged temperature and entropy more inline with that of run RT000 V000 P0, offsetting these curves such that they lie above those of the other inclination angles. Changes in the volume averaged temperature and entropy between the other wind angles are found to be modest in comparison. Additionally, changes in the volume averaged cocoon temperature and entropy (alongside the pressure and density) as a result of changes in the impact parameter are minimal (compare the solid red curve and dash-dotted blue and magenta curves in Fig. 14). The cocoon parameters therefore depend more strongly on the nature of the wind (velocity and inclination angle) rather than the impact parameter between the jets -at least for the values of P simulated.

0 E N S T RO P H Y
In order to investigate the turbulence associated with the interaction of the jets, we look at the enstrophy density (r), which is given by the mean of the magnitude of the vorticity squared: The vorticity at position r, ω(r), is given by ∇ × v(r) and can be used to quantify the extent of the rotation in a turbulent fluid flow. We note that the turbulence measured here is highly dependent on grid resolution as we do not include a sub-grid prescription for turbulence in these simulations. The refinement criteria used in the simulations are based on the second derivative of the density, ensuring that regions with high-density contrasts, such as those in the vicinity of shocks, are maximally refined. Furthermore, maximum refinement is enforced in regions that contain high values of the jet tracer fractions, ensuring that the bulk jet flow is highly refined. (Maximum refinement is always enforced at the location of the jet nozzles to ensure that the correct momentum flux of the jet is injected at each time-step.) As the jets propagate thorough the ambient medium, they slow as a result of their bulk kinetic energy being dissipated into turbulence and shock waves. The enstrophy can be used as a measure of turbulence, as has been done to investigate the turbulent nature of the ICM in galaxy clusters by Porter, Jones & Ryu (2015), Wittor et al. (2017), and Vazza et al. (2017). Discussions concerning enstrophy and the generation of enstrophy, with regard to AGN jets, however, remain comparatively sparse. We note that our first paper (Molnar et al. 2017) appears to be the first publication, as far as we have found, in which the enstrophy is used to discuss the turbulence in simulated AGN jets. Following on from our first paper, we conduct a similar analysis of the enstrophy and its generation, expanding the discussion in Molnar et al. (2017) to include the effects of the ICM wind on the colliding jets. Fig. 15 shows the enstrophy projected along the z-axis for a selection of the simulations. The figure shows that the regions with highest enstrophy, and thus the regions with the largest rates of energy dissipation, are close to the jet axis. The more powerful bottom two-sided jet in each of the simulations displays a higher enstrophy along its length compared to the weaker jet. The enstrophy of both the top and bottom jets is enhanced at distances prior to, and immediately following, the deflection point for regions of the jets directly exposed to the wind. The largest source of enstrophy in the simulations typically corresponds to the merging jets (the jets propagating to the right of the simulation domain), and higher values of the enstrophy are present in the merged jet when the wind is aligned such that it is almost head-on to the direction of propagation of the merging jets, as in RT225 V120 P0. The enhancement in enstrophy as a result of the direct contact between the jet and the inner cocoon edge is also significant. Fig. 16 shows the total enstrophy in the cocoon region (which includes the bulk jet flow in its definition) as a function of evolution time for runs RT000 V000 P0 (solid black curve), RT130 V120 P0 (dotted red curve), and RT130 V120 P1 (dashed-dot blue curve). The enstrophy shown is normalized by its value in the cocoon region at t = 0.5 Myr (this corresponds to the second plotfile in each run). For t < 4 Myr, the enstrophy evolves similarly in for each run. Beyond t = 4 Myr, the curves begin to diverge and the enstrophy grows more slowly for run RT130 V120 P0 (conducted with impact parameter P = 0 and wind velocity of 1200 km s −1 ) for the majority of the remaining evolution. Beyond t = 27 Myr, the cocoon enstrophy is largest for run RT000 V000 P0 conducted without the ICM wind.

Enstrophy generation
In order to explore how enstrophy is generated in the simulations and identify the motions that comprise the main sources of enstrophy, we require an expression for the temporal evolution of the enstrophy. The temporal evolution of the enstrophy is given by   Total enstrophy in the cocoon as a function of evolution time for simulations RT000 V000 P0 (no wind; black curve), RT130 V120 P0 (wind velocity of 1200 km s −1 and wind angle of 130 • ; dotted red curve), and RT130 V120 P1 (dash-dotted blue curve; impact parameter P = 1, wind parameters as in RT130 V120 P0). The enstrophy is normalized by its value in the second plotfile (t = 0.5 Myr) for each simulation and is largest for the simulation conducted without the wind from t = 22 Myr. F comp , F baro , F str , and F diss describe the sources and sinks of the enstrophy and govern the evolution of the enstrophy. As per the discussion in Wittor et al. 2017, these terms can be outlined as follows.
F advec describes the enstrophy generation due to advective motions and is given by The term describes the conservative advection of enstrophy across the simulation domain. The compressive term F comp accounts for the amplification of enstrophy as a result of shocks, alongside reversible compressions and rarefractions (Wittor et al. 2017) and is given by The baroclinic term F baro is a result of the misalignment between the gradients of the gas density and pressure and generates enstrophy in baroclinic flows (Porter et al. 2015) -flows in which the gas pressure is not solely a function of density. The baroclinic term is given by The term for the stretching motions, F str , accounts for enstrophy generation via vortex stretching (Wittor et al. 2017) and is given by Vortex stretching dominates the enstrophy generation in cluster interiors, although the other terms also contribute substantially ). The dissipation term F diss accounts for the viscous dissipation of the flow and is given by where ν is the kinematic viscosity and G = (1/ρ)∇ρ · S, where S is the strain tensor (see Mee & Brandenburg 2006). Since we do not explicitly introduce viscosity in these simulations, the kinematic viscosity ν is formally 0. However, as described by Porter et al. (2015), the numerical dissipation inherent in the simulations can mimic viscosity on length scales comparable to the grid resolution. Quantifying this 'numerical viscosity' due to the numerical dissipation in the simulations is non-trivial, and so we neglect the dissipative term F diss when computing the evolution of the enstrophy.
We compute the enstrophy growth terms for the cocoon region. Fig. 17 shows the enstrophy growth terms for simulations RT000 V000 P0 (solid black curves), RT130 V120 P0 (dotted red curves), and RT130 V120 P1 (dashed-dotted blue curves) integrated over the cocoon volume. As we are interested in the relative differences between terms, all curves are normalized by the same arbitrary factor. The stretching term F str (bottom righthand panel of Fig. 17) dominates the other growth terms during the evolution for each simulation, and the shape of the curve most closely resembles that of the total enstrophy shown in Fig. 16 for each run. This implies that vortex stretching is the driving mechanism for enstrophy generation in the simulations. The compressive term F comp (Fig. 17, bottom left-hand panel) and thus shocks in the simulations increasingly contribute to the enstrophy generation as the evolution time increases. The baroclinic and advective terms, F baro and F advec (Fig. 17, top left-hand and right-hand panels, respectively), contribute less significantly to the enstrophy generation. Vortex stretching and baroclinic motions contribute more to the enstrophy generation in Run RT000 V000 P0 conducted Figure 17. Enstrophy growth terms F baro (top left) and F advec (top right), F comp (bottom left) and F str (bottom right) for runs RT000 V000 P0 (solid black curve), RT130 V120 P0 (dotted red curve), and RT130 V120 P1 (dash-dotted blue curve). See Section 10.1 for more information. The curves are normalized by the same arbitrary factor. The stretching term F str dominates the enstrophy growth, and there is a significant contribution to enstrophy growth from the compression term F comp . Contributions to enstrophy via advection (F advec ) and baroclinic motions (F baro ) are comparatively negligible. without the ICM wind, which consistently displays the largest F str and F baro from evolution times beyond t ≈ 2.3 Myr. Compressive motions become more significant for runs RT130 V120 P0 and RT130 V120 P1 compared to run RT000 V000 P0 without the wind from t ≈ 11.5 Myr. Whilst the least significant term F advec remains similar for all three runs, indicating that advective motions are least significant for enstrophy generation.

C O N C L U S I O N S A N D F U T U R E W O R K
The work presented here comprises the first hydrodynamic simulations of two bipolar jets propagating in an external wind for application to the binary AGN 3C 75. We aimed to investigate the effect of wind direction, velocity, and jet impact parameter on the evolution, morphology, and energetics of the jets in order to begin constraining the interactions/conditions behind the striking morphology of 3C 75. This work extends the numerical study presented in Molnar et al. (2017) where we investigate the collisions between jets of different Mach numbers and collisional impact parameters in a static, uniform ambient medium. The main conclusions of the work presented here are as follows: (i) A wind velocity of 1200 km s −1 is sufficient to deflect both bipolar jets into the direction of the wind. Simulations with wind angles θ w ∼ 100 • −150 • and a wind velocity v w = 1200 km s −1 provide the closest match to the morphology observed in 3C 75 for the wind parameters simulated. For these wind angles, both jets propagate within a common cocoon on both sides of the source, the resulting morphology of which is highly perturbed by the wind.
(ii) Different wind orientations, wind velocity, and collisional jet impact parameters generate different shock features, pressure, and density variations within the shocked shell of ICM surrounding the jets. These variations will lead to different regions of enhanced X-ray surface brightness within the shocked shell. Simulations in which cross winds are simulated predict regions of high density and pressure in the inner regions of the shocked shell. These density enhancements arise as the bulk jet flow is pushed against the inner edge of the cocoon by the cross wind. If variations in the X-ray surface brightness of the shocked shell region are resolved in highexposure X-ray observations, they may be able to shed light on the nature of the interaction between the jets.
(iii) The collision between the jets results in a standing shock that remains stationary throughout the evolution, regardless of the wind parameters. The location of the standing shock corresponds to a local region of high pressure and density and could be a site of strong localized particle acceleration and enhanced synchrotron emission. If entrained thermal material is advected with the relativistic jet plasma, then the collision region may be associated with an increase in the X-ray surface brightness.
(iv) The wind and the collision can cause the jets to disrupt into kiloparsec scale filaments. This would lead to neighbouring regions of high-and low-surface brightness radio emission, or regions of Faraday rotation banding. If the filaments oscillate and cross each other, then they may produce a helical topology on kiloparsec scales, such as that observed in 3C 75.
(v) If the impact parameter is sufficiently large that the jets 'miss' one another (P = 1 and 3), then the wind is not sufficient to merge their trajectories on the right of the source and they appear to cross each other in projection. Furthermore, when the bipolar jets are simulated separately (runs RT130 V120 top and RT130 V120 bot), their location of maximum curvature does not happen at the same location in the atmosphere. Both the results from the impact parameter runs and single jet simulations therefore suggest that some direct contact between the bulk jet flows is required, in order for the jets to bend in tandem on the northwestern side of 3C 75. Whether the collision is head-on (P = 0) or corresponds to a small finite impact parameter (0 < P < 1) remains to be determined.
(vi) The wind acts to increase the internal, kinetic, and potential energy in the cocoon region at a given evolution time. Changing the wind velocity can induce significant changes in cocoon energetics, with faster wind velocities resulting in higher cocoon energies for all energy components. Changes in the energetics of the cocoon as a result of varying the wind direction are less significant, while changes in the cocoon energetics due to the impact parameter of the colliding jets are comparatively negligible.
(vii) The presence of the wind increases the level of entrainment in the cocoons. The level of entrainment shows a strong dependence on wind velocity and orientation. Cross winds result in higher levels of entrainment as prolonged contact between the bulk flow and inner edge of the cocoon can further disrupt the contact discontinuity, promoting the development of disruptive Kelvin-Helmholtz instabilities at the contact discontinuity. Non-zero impact parameters also lead to a significant increase in the level of entrainment due to the prolonged contact between the bulk jet flow of the weaker top jet with the inner edge of its cocoon. The amount of cold, dense, low entropy shocked ICM entrained into the cocoon can also influence the global parameters of the cocoons -simulations conducted without the wind entrain the least ICM and produce cocoons that are less dense, hotter, and contain material of higher entropy compared to the simulations with the wind.
(viii) Enstrophy is enhanced in the simulation conducted without the wind. Vortex stretching is the driving mechanism for enstrophy generation in the cocoons (which includes the bulk jet flow in its definition), while the contribution from compressive motions is also significant. Baroclinic and advective motions are less significant in comparison.
In this work, we simulate the jets in a purely hydrodynamic framework. The inclusion of a dynamically important magnetic field may affect the growth of instabilities and thus the advanced speed of the jet into the ambient medium. This could have a significant impact on both the degree of curvature and the stability of the jet (and the development of filamentary structures) in the presence of a wind. The effect of magnetic fields on the collision and deflection of the jets is beyond the scope of this work, and we leave this for future investigation.
Throughout this work, we have also assumed that both bipolar jets have the same density, an uneven distribution of the total kinetic power of the source between the two bipolar jets and that the jets 'switch-on' simultaneously. All these assumptions will have implications on the resulting morphology of the jets, and we leave this investigation for future work.
In order to conduct a more in-depth simulation study relevant to 3C 75, we need better constraints on the parameters of the ICM. Chandra observations with an exposure time much greater than those presented in Hudson et al. (2006) are required in order to obtain a detailed description of density and temperature profiles for the cluster environment in Abell 400. Further VLA radio observations will also allow better constraints for some of the source properties, such as the projection angle of 3C 75. Polarization studies of 3C 75 may confirm the presence of a jet-on-jet collision within the source.

AC K N OW L E D G E M E N T S
This work was carried out using the computational facilities of the Advanced Computing Research Centre, University of Bristol. 4 The FLASH (M)HD code used in this work was in part developed by the DOE NNSA-ASC OASCR Flash Center at the University of Chicago. We thank Josie Rawes for generating the VLA radio map shown in Fig. 1. G. M. acknowledges support from the Netherlands Organisation for Scientific Research (NWO) VICI award, grant no. 639.043.513 and a Netherlands Research School for Astronomy (NOVA), Virtual Institute of Accretion (VIA) postdoctoral fellowship. This research has made use of the VisIt visualization package (Childs et al. 2012) and the YT analysis toolkit (Turk et al. 2011).
Processed data for Figs 11-14 and Figs 16-17 presented in this work will be made available at doi:10.5281/zenodo.3747409.