This white paper describes the science case for Very Long Baseline Interferometry (VLBI) and provides suggestions towards upgrade paths for the European VLBI Network (EVN). The EVN is a distributed long-baseline radio interferometric array, that operates at the very forefront of astronomical research. Recent results, together with the new science possibilities outlined in this vision document, demonstrate the EVN's potential to generate new and exciting results that will transform our view of the cosmos. Together with e-MERLIN, the EVN provides a range of baseline lengths that permit unique studies of faint radio sources to be made over a wide range of spatial scales. The science cases are reviewed in six chapters that cover the following broad areas: cosmology, galaxy formation and evolution, innermost regions of active galactic nuclei, explosive phenomena and transients, stars and stellar masers in the Milky Way, celestial reference frames and space applications. The document concludes with identifying the synergies with other radio, as well as multi-band/multi-messenger instruments, and provide the recommendations for future improvements. The appendices briefly describe other radio VLBI arrays, the technological framework for EVN developments, and a selection of spectral lines of astrophysical interest below 100 GHz. The document includes a glossary for non-specialists, and a list of acronyms at the end.

Hyper-accretion occurs when the gas inflow rate onto a black hole (BH) is so high that the radiative feedback cannot reverse the accretion flow. This extreme process is a promising mechanism for the rapid growth of seed BHs in the early universe, which can explain high-redshift quasars powered by billion solar mass BHs. In theoretical models, spherical symmetry is commonly adopted for hyper-accretion flow, however its sustainability on timescales corresponding to the BH growth has not been addressed yet. Here we show that stochastic interactions between the ionizing radiation from the BH and non-uniform accretion flow can lead to a formation of rotating gas disk around the BH. Once the disk is formed, the supply of gas to the BH preferentially occurs via biconical dominated accretion flow (BDAF) perpendicular to the disk, avoiding the centrifugal barrier of the disk. BDAFs from opposite directions collide in the vicinity of the BH supplying high-density, low angular momentum gas to the BH, whereas most of the gas with non-negligible angular momentum is deflected to the rotationally supported outflowing decretion disk. The disk becomes reinforced progressively as more mass from the biconical flow transfers to the disk and some of the outflowing gas from the disk is redirected to the biconical accretion funnels through a meridional structure. This axisymmetric hydrodynamic structure of a BDAF and decretion disk continues to provide uninterrupted flow of high-density gas to the BH.

The merging of supermassive black holes (SMBHs) is a direct consequence of our hierarchical picture of galaxy evolution. It is difficult to track the merging process of SMBHs during mergers of galaxies as SMBHs are naturally difficult to observe. We want to characterise and confirm the presence of two independent active galactic nuclei (AGN) separated by a few kiloparsec in seven strongly interacting galaxies previously selected from the Sloan Digital Sky Survey (SDSS) as Seyfert-Seyfert pairs based on emission-line ratio diagnostics. Optical slit spectra taken with MODS at the Large Binocular Telescope (LBT) are presented to infer the detailed spatial distribution of optical emission lines, and their line ratios and AGN signatures with respect to the host galaxies, thereby quantifying the impact of beam smearing and large fibre apertures on the spectra captured by the SDSS. We find that at most two of the seven targets actually retain a Seyfert-Seyfert dual AGN, whereas the others may be more likely powered by post-AGB stars in retired galaxies or through shocks in the ISM based on spatially resolved optical line diagnostics. The major cause of this discrepancy is a bias caused by the spillover of flux from the primary source in the secondary SDSS fibre which can be more than an order of magnitude at <3" separations. Previously reported extremely low X-ray-to-[\ion{O}{iii}] luminosity ratios may be explained by this misclassification, as can heavily obscured AGN for the primaries. We also find that the nuclei with younger stellar ages host the primary AGN. Studies of close dual AGN selected solely from fibre-based spectroscopy can create severe biases in the sample selection and interpretation of the results. Spatially resolved spectroscopy should ideally be used in the future to characterise such compact systems together with multi-wavelength follow-up observations.

Dynamical friction is thought to be a principal mechanism responsible for orbital evolution of massive black holes (MBHs) in the aftermath of galactic mergers and an important channel for formation of gravitationally bound MBH binaries. We use 2D radiative hydrodynamic simulations to investigate the efficiency of dynamical friction in the presence of radiative feedback from an MBH moving through a uniform density gas. We find that ionizing radiation that emerges from the innermost parts of the MBH's accretion flow strongly affects the dynamical friction wake and renders dynamical friction inefficient for a range of physical scenarios. MBHs in this regime tend to experience positive net acceleration, meaning that they speed up, contrary to the expectations for gaseous dynamical friction in absence of radiative feedback. The magnitude of this acceleration is however negligibly small and should not significantly alter the velocity of MBHs over relevant physical timescales. Our results suggest that suppression of dynamical friction is more severe at the lower mass end of the MBH spectrum which, compounded with inefficiency of the gas drag for lower mass objects in general, implies that $< 10^7$ solar mass MBHs have fewer means to reach the centers of merged galaxies. These findings provide formulation for a sub-resolution model of dynamical friction in presence of MBH radiative feedback that can be easily implemented in large scale simulations.

The quest for binary and dual supermassive black holes (SMBHs) at the dawn of the multi-messenger era is compelling. Detecting dual active galactic nuclei (AGN) -- active SMBHs at projected separations larger than several parsecs -- and binary AGN -- probing the scale where SMBHs are bound in a Keplerian binary -- is an observational challenge. The study of AGN pairs (either dual or binary) also represents an overarching theoretical problem in cosmology and astrophysics. The AGN triggering calls for detailed knowledge of the hydrodynamical conditions of gas in the imminent surroundings of the SMBHs and, at the same time, their duality calls for detailed knowledge on how galaxies assemble through major and minor mergers and grow fed by matter along the filaments of the cosmic web. This review describes the techniques used across the electromagnetic spectrum to detect dual and binary AGN candidates and proposes new avenues for their search. The current observational status is compared with the state-of-the-art numerical simulations and models for formation of dual and binary AGN. Binary SMBHs are among the loudest sources of gravitational waves (GWs) in the Universe. The search for a background of GWs at nHz frequencies from inspiralling SMBHs at low redshifts, and the direct detection of signals from their coalescence by the Laser Interferometer Space Antenna in the next decade, make this a theme of major interest for multi-messenger astrophysics. This review discusses the future facilities and observational strategies that are likely to significantly advance this fascinating field.

Accretion of gas and interaction of matter and radiation are at the heart of many questions pertaining to black hole (BH) growth and coevolution of massive BHs and their host galaxies. To answer them it is critical to quantify how the ionizing radiation that emanates from the innermost regions of the BH accretion flow couples to the surrounding medium and how it regulates the BH fueling. In this work we use high resolution 3-dimensional (3D) radiation-hydrodynamic simulations with the code Enzo, equipped with adaptive ray tracing module Moray, to investigate radiation-regulated BH accretion of cold gas. Our simulations reproduce findings from an earlier generation of 1D/2D simulations: the accretion powered UV and X-ray radiation forms a highly ionized bubble, which leads to suppression of BH accretion rate characterized by quasi-periodic outbursts. A new feature revealed by the 3D simulations is the highly turbulent nature of the gas flow in vicinity of the ionization front. During quiescent periods between accretion outbursts, the ionized bubble shrinks in size and the gas density that precedes the ionization front increases. Consequently, the 3D simulations show oscillations in the accretion rate of only ~2-3 orders of magnitude, significantly smaller than 1D/2D models. We calculate the energy budget of the gas flow and find that turbulence is the main contributor to the kinetic energy of the gas but corresponds to less than 10% of its thermal energy and thus does not contribute significantly to the pressure support of the gas.

We examine radiation-regulated accretion onto intermediate-mass and massive black holes (BHs) embedded in a bulge component. Using spherically symmetric one-dimensional radiation-hydrodynamics simulations, we track the growth of BHs accreting from a cold, neutral gas reservoir with temperature T=10^4 K. We find that the accretion rate of BHs embedded in bulges is proportional to r_{B,eff}/r_B, where r_{B,eff} is the increased effective Bondi radius that includes the gravitational potential of the bulge, and r_B is the Bondi radius of the BH. The radiative feedback from the BH suppresses the cold accretion rate to ~1 percent of the Bondi rate when a bulge is not considered. However, we find that the BH fueling rate increases rapidly when the bulge mass M_bulge is greater than the critical value of 10^6 M_sun and is proportional to M_bulge. Since the critical bulge mass is independent of the central BH mass M_{BH}, the growth rate of BHs with masses of 10^2, 10^4, and 10^6 M_sun exhibits distinct dependencies on the bulge-to-BH mass ratio. Our results imply that light seed BHs (<= 10^2 M_sun) which might be the remnants of the Pop III stars, cannot grow through accretion coevally with the early assembly of the bulge of the host galaxies until the bulge reaches the critical mass. However, massive BH seeds (>= 10^5 M_sun) that may form via direct collapse, are more likely to be embedded in a supercritical bulge and thus can grow efficiently coupling to the host galaxies and driving the early evolution of the M_{BH}-$\sigma$ relationship.

We used a model of a relativistic accretion disk around a supermassive black hole (SMBH), based on ray-tracing method in the Kerr metric, to study the variations of the composite Fe K$\alpha$ line emitted from two accretion disks around SMBHs in a binary system. We assumed that the orbit of such a binary is approximately Keplerian, and simulated the composite line shapes for different orbital elements, accretion disk parameters and mass ratios of the components. The obtained results show that, if observed in the spectra of some SMBH binaries during their different orbital phases, such composite Fe K$\alpha$ line profiles could be used to constrain the orbits and several properties of such SMBH binaries.

The role played by magnetic fields in the intracluster medium (ICM) of galaxy clusters is complex. The weakly collisional nature of the ICM leads to thermal conduction that is channelled along field lines. This anisotropic heat conduction profoundly changes the stability of the ICM atmosphere, with convective stabilities being driven by temperature gradients of either sign. Here, we employ the Athena magnetohydrodynamic code to investigate the local non-linear behavior of the heat-flux driven buoyancy instability (HBI), relevant in the cores of cooling-core clusters where the temperature increases with radius. We study a grid of 2-d simulations that span a large range of initial magnetic field strengths and numerical resolutions. For very weak initial fields, we recover the previously known result that the HBI wraps the field in the horizontal direction thereby shutting off the heat flux. However, we find that simulations which begin with intermediate initial field strengths have a qualitatively different behavior, forming HBI-stable filaments that resist field-line wrapping and enable sustained vertical conductive heat flux at a level of 10--25% of the Spitzer value. While astrophysical conclusions regarding the role of conduction in cooling cores require detailed global models, our local study proves that systems dominated by HBI do not necessarily quench the conductive heat flux.