Magnetohydrodynamic modelling of star-planet interaction and associated auroral radio emission

We present calculations of auroral radio powers of magnetised hot Jupiters orbiting Sun-like stars, computed using global magnetohydrodynamic (MHD) modelling of the magnetospheric and ionospheric convection arising from the interaction between the magnetosphere and the stellar wind. Exoplanetary auroral radio powers are traditionally estimated using empirical or analytically-derived relations, such as the Radiometric Bode's Law (RBL), which relates radio power to the magnetic or kinetic energy dissipated in the stellar wind-planet interaction. Such methods risk an oversimplification of the magnetospheric electrodynamics giving rise to radio emission. As the next step toward a self-consistent picture, we model the stellar wind-magnetosphere-ionosphere coupling currents using a 3D MHD model. We compute electron-cyclotron maser instability-driven emission from the calculated ionospheric field-aligned current density. We show that the auroral radio power is highly sensitive to interplanetary magnetic field (IMF) strength, and that the emission is saturated for plausible hot Jupiter Pedersen conductances, indicating that radio power may be largely independent of ionospheric conductance. We estimate peak radio powers of $10^{14}$ W from a planet exposed to an IMF strength of $10^3$ nT, implying flux densities at a distance of 15 pc from Earth potentially detectable with current and future radio telescopes. We also find a relation between radio power and planetary orbital distance that is broadly consistent with results from previous analytic models of magnetosphere-ionosphere coupling at hot Jupiters, and indicates that the RBL likely overestimates the radio powers by up to two orders of magnitude in the hot Jupiter regime


INTRODUCTION
The direct detection of exoplanets in large parts of the electromagnetic spectrum is hindered by the high luminosity contrast ratio between the star and planet. Evidence from the Solar system planets, however, indicates that the radio waveband presents a luminosity ratio much more conducive to direct detection, with non-thermal emission from Jupiter of similar intensity to Solar radio bursts (Zarka 2007). Historically, the search for exoplanetary radio emission has focussed primarily on Jupiter-like exoplanets in close orbit (3 -10 stellar radii) around their parent star. Such planets, given the moniker 'hot Jupiters', have been targeted by several studies examining the possibility of detectable auroral E-mail: jdn4@leicester.ac.uk (KTS) radio emission (e.g. Farrell et al. 1999;Zarka et al. 2001;Zarka 2007;Grießmeier et al. 2004;Griessmeier et al. 2005Griessmeier et al. , 2007Lazio et al. 2004;Jardine & Cameron 2008;Nichols 2011Nichols , 2012Hallinan et al. 2012;Vidotto & Donati 2017). Although the search for radio emission has not been limited exclusively to hot Jupiters, the plausibly strong magnetic fields and intense stellar wind conditions present at these objects is thought to be favourable to the generation of radio emission through star-planet interaction. Auroral radio emission from exoplanets is envisaged to be generated by the electron-cyclotron maser instability (ECMI), the same mechanism responsible for driving radio emission from auroras at magnetized Solar system planets (Wu & Lee 1979;Treumann 2006;Imai et al. 2008;Lamy et al. 2011). This intense, coherent electromagnetic radiation is emitted at the local cyclotron frequency, and therefore Jupiter-like exoplan-ets are the prime candidates for directly detectable emission, since their potentially high intrinsic magnetic field strengths (∼ B Jup ) are required to produce emission above the Earth's ionospheric cutoff frequency of ∼ 10 MHz.
Many previous studies estimating expected exoplanetary radio emission employ the Radiometric Bode's Law (RBL), an empirical scaling relation extrapolated from Solar system measurements between incident Poynting or kinetic energy flux and emitted radio power (Farrell et al. 1999;Zarka et al. 2001;Zarka 2007;Lazio et al. 2004). Estimates of the radio power from hot Jupiter auroras based on the RBL range between 10 14 -10 16 W, up to five orders of magnitude stronger than Jupiter's equivalent ECMIdriven emission (Zarka 2007), implying radio fluxes which may be detectable with the next generation of radio telescopes. Two primary processes are assumed to mediate radio emission driven by star-planet interaction: either Alfvén waves, such as mediates emission in the sub-Aflvénic Io-Jupiter interaction; or magnetic reconnection, such as occurs at Earth's magnetosphere. In considering radio emission propagated by Alfvén waves, Saur et al. (2013) and Turnpenney et al. (2018) estimated total radiated energy fluxes from exoplanets of up to 10 19 W, translating to radio powers of 10 17 W assuming an ECMI efficiency of 1 per cent. Nichols & Milan (2016), using an analytic model, considered an Earth-type Dungey cycle process of magnetic reconnection, computing ionospheric field-aligned currents (FACs) and resulting radio powers arising from ionospheric convection. They found that saturation of the convectioninduced electric potential limited the dissipated power, and predicted auroral radio emission from hot Jupiters ∼2 orders of magnitude smaller than RBL-based predictions. In this paper we use a global 3D magnetohydrodynamic (MHD) model to calculate the magnetosphere-ionosphere coupling currents arising from the interaction between Sun-like stars and hot Jupiters, and hence estimate the resulting auroral radio emission. By using a numerical global MHD model, this study extends the analytic work of Nichols & Milan (2016) to enable a self-consistent calculation of the FACs at exoplanets. The numerical model, as used in this study, takes a set of input boundary conditions and generates a 2D map of the ionospheric FAC density distribution, from which radio power is calculated in post-processing. The overall power is expected to be strongly influenced externally, by the interplanetary magnetic field (IMF) strength, and internally, by the ionospheric Pedersen conductance. Therefore, we run simulations across a broad range of these two parameters, comparing the results with those of Nichols & Milan (2016). This paper begins with an overview of the theoretical background relevant to the MHD model, along with the formulation used to compute radio powers from the fieldaligned current output of the model. Then follows a description of the method employed, before results of the modelling work are presented and discussed. A case intermediate between Earth and hot Jupiter exoplanets is first studied, before cases more appropriate to exoplanets are examined.

Magnetohydrodynamic model
Global 3D MHD simulations, based on first principle physics, are a powerful tool for modelling the dynamics and evolution of magnetic fields and plasma flows in space weather and astrophysical phenomena. In use since the 1980s, early models relied on techniques such as finitedifference methods to solve a system of discretized MHD equations (Van Leer 1979;Leboeuf et al. 1981;Wu et al. 1981;Brecht et al. 1982). Computational solutions of MHD equations require discretization of the MHD equations, which inherently introduces errors into the solution. Modern MHD models use advanced solution techniques which improve the efficiency of the solution and minimise such discretization errors . Such methods rely on approximations to solve a Riemann problem, the form of initial value problem presented in MHD numerical analysis over a finite volume. The 3D MHD solver used in this work is the 'Block Adaptive Tree Solarwind Roe Upwind Scheme' (BATS-R-US) software first outlined by Powell et al. (1999), and developed at the University of Michigan. The computational scheme of BATS-R-US is based on the same elements used in many state-of-the-art MHD codes, and this section describes each of those elements in the scheme.
A governing set of 3D ideal MHD equations is first defined. Various forms of these equations are expressible, and the form chosen is dictated by factors which will ultimately aid in the computational solution of these equations. The set of ideal MHD equations used in BATS-R-US are expressed in a gasdynamic conservative form ) where I is the identity matrix and the total gasdynamic energy E gd , is given by where γ is the ratio of the specific heats. This equation set contains an expression for the conservation of mass (equation 1), conservation of momentum (equation 2), an expression of Faraday's law (equation 3), and an energy equation (4). These partial differential equations are manipulated into a non-dimensional, symmetrizable form, the full details of which can be found in Powell et al. (1999). A computational domain is divided into cells over which the MHD equations are integrated. From the symmetrizable form Powell (1994) showed that it is possible to derive a Roe-type approximate Riemann solver for the 3D equations. First described by Roe (1981), this is a method for solving partial differential equations by estimating the flux at the interface between two computation cells in some finite-volume discretized domain. Such solvers are required in magnetohydrodynamics, since iterative solutions are costly, and therefore approximations must be made. The computation domain is divided into a grid of Cartesian cells, and the cells are structured into blocks typically consisting of between 4 × 4 × 4 and 12 × 12 × 12 individual cells. The block-adaptive technology of BATS-R-US allows the computational grid to be adapted based on prespecified physical criteria, such that blocks can be refined in regions where interesting physical features emerge. Adaptive mesh refinement is extremely effective when the problem being treated contains disparate length scales, and also removes any initial grid-based bias in the solution.
The Space Weather Modelling Framework (SWMF) is a software package which integrates several different physics domains extending from the solar surface to the planetary upper altmosphere (Tóth et al. 2005). BATS-R-US is used in a number of these components where MHD solutions are required, and physically meaningful combinations of components can be coupled together to study a wide variety of space weather events and phenomena (Tóth et al. 2012). While developed originally to study Sun-Earth events, the SWMF has since been adapted and applied to other solar system planets, satellites and comets (e.g. Jia et al. 2019;Jia & Kivelson 2016;Tóth et al. 2016;Huang et al. 2016), and may reasonably be used for the study of extrasolar astrophysical systems where the physics domains are appropriate, with some adaptation where may be required.

Field-aligned current and radio power
This study utilizes the Global Magnetosphere (GM) component of the SWMF, coupled with the Ionospheric Electrodynamics (IE) component (Ridley et al. 2004). The GM domain constructs the magnetic environment and plasma dynamics around the planet, and contains features such as the bow shock, magnetopause, and magnetotail. Upstream boundary condition for the GM component can be obtained from coupling with the Inner Heliosphere component of the SWMF, but in this work are simply input into the model based on reasonable values, as will be discussed below. Currents from the GM component are mapped down along magnetic field lines to provide the field-aligned current boundary conditions for the IE component. The domain of the IE component is a height-integrated spherical surface. Formally, this component is a two-dimensional electric potential solver, which computes conductances and particle precipitation from FACs. The process can be summarized as follows: 1) Field-aligned currents are calculated by ∇ × B at 3.5 R P , a value also employed by Ridley et al. (2004), where B is the local magnetic field; 2) The currents are then mapped down along field lines to a nominal ionospheric altitude of ∼ 110 km using the planetary dipolar field, and are scaled by the ratio B I /B 3.5 , where B I and B 3.5 are the magnetic field strengths at the ionosphere and 3.5 planetary radii respectively. 3) Next, a height-integrated ionospheric conductance map is generated and the electric potential is calculated, which is then mapped out along magnetic field lines to the simulation's inner boundary at 2.5 R P , where flow velocities and electric fields are prescribed (Ridley et al. 2001(Ridley et al. , 2004. Of the several variables output from the IE component, this work is principally concerned with the ionospheric FAC density j | | , and the cross-polar cap potential (CPCP) Φ CPC outputs. By integrating the total upward or downward FAC density output from the IE component over one hemisphere, the total current, I tot , flowing into or out of the ionosphere is obtained by where θ and φ are the conventional spherical coordinates of colatitude and azimuth respectively. Total auroral power is also calculated in post-processing by integrating precipitating electron energy flux E f over one hemisphere: Qualitatively, the precipitating electron energy flux is the kinetic energy carried by the downward-flowing electrons associated with the upward field-aligned current. The maximum FAC which can be carried by unaccelerated electrons in an isotropic Maxwellian velocity space distribution is given by where e and m e are the electron charge and mass respectively, and W th and n are are the thermal energy and number density of the magnetospheric electron source population respectively, for which we employ canonical jovian values established from Voyager measurements of W th = 2.5 keV and n = 0.01 cm −3 throughout this work (e.g. Scudder et al. 1981). We note that these values may differ significantly at hot Jupiters, although no reliable estimates exist at present. The source plasma number density is in relation to the evacuated auroral field lines, and is therefore expected to be much reduced from the ambient plasma density. Qualitatively, the effect of varying these parameters, however, will be to increase the precipitating electron energy flux and auroral power with increasing plasma thermal energy, and with decreasing plasma density. The jovian values are employed here in the absence of any information for realistic values at hot Jupiters, although future work may investigate a range of these parameter values to determine quantitatively the effect on auroral radio emission. In general, the FACs will be larger than can be carried solely by unaccelerated electrons, and therefore must be driven by a field-aligned electric potential. In common with previous works computing intense exoplanetary and ultracool dwarf auroral radio emissions (Nichols 2011(Nichols , 2012Nichols et al. 2012;Turnpenney et al. 2017) , we employ Cowley's (2006) relativistic extension of Knight's (1973) current-voltage relation for parallel electric fields given by where Φ min is the minimum field-aligned voltage required to drive the current j i in the ionosphere, and c is the speed of light in a vacuum. The corresponding precipitating electron energy flux is given by where E 0 is the maximum unaccelerated electron energy flux, corresponding to equation (8), given by Assuming, in common with observations of Jupiter and Saturn, an ECMI efficiency of 1 per cent (Gustin et al. 2004;Clarke et al. 2009;Lamy et al. 2011), the emitted auroral radio power is Finally, the spectral flux density is calculated by where s is the distance to the exoplanetary system from the Earth, ∆ν is the radio emission bandwidth, and a beaming angle of 1.6 sr is assumed on the basis of Jupiter's observed ECMI emission beaming angle (Zarka et al. 2004).
Since cyclotron maser emission is generated at the local electron-cyclotron frequency, the bandwidth is the difference between the magnetic field strength at the location of the field-aligned potential and the ionosphere. This formulation assumes that the potential is located high enough up the field line from the ionosphere that the field strength there is much smaller than the field strength in the ionosphere. This assumption is valid beyond a few planetary radii owing to the r −3 dependence for a dipole planetary field. We therefore assume that the bandwidth ∆ν is equal to the electron cyclotron frequency in the polar ionosphere and hence given by where B i is the ionospheric magnetic field strength.

METHODOLOGY
Input parameters were chosen for the SWMF to simulate a magnetised hot Jupiter interacting with the IMF and stellar wind of a solar-type star. A planet of Jupiter mass (1.9 × 10 27 kg), radius (71, 492 km), and equatorial magnetic field strength (426, 400 nT) was specified, with a dipole magnetic field aligned with the planetary spin axis. The orientation of the magnetic field is opposite to that of Jupiter, i.e. the planetary field is pointing northward at the equator. Each simulation was run with a planetary plasma density of 10 7 cm −3 , and a temperature of 8, 000 K at the inner boundary of the simulation, values consistent with those based on modelling of atomic hydrogen, heavy atoms and ions surrounding hot Jupiters (Muñoz 2007;Koskinen et al. 2013;Shaikhislamov et al. 2016). The incident plasma velocity of the simulations was set to 250 km s −1 . This value represents the impinging plasma flow velocity taking into account both the stellar wind outflow and the Keplerian orbital velocity of the planet in the hot Jupiter regime. We note that although in reality the incident plasma velocity will be predominantly azimuthal due to the high orbital velocity of planet, in this work this incident plasma velocity is prescribed as an input in the model, such that in the results that follow (i.e. Figs 1 and 2) the X-axis is aligned with the incident plasma velocity. This transposition is purely for convenience of modelling, and does not affect the resulting auroral radio power calculations.
Since the stellar wind plasma temperature close to the star in the orbital distance regime of hot Jupiters is approximately the coronal value, a stellar wind plasma temperature input of 2 MK was employed throughout this study. The stellar wind density was set to 10 4 cm −3 in all simulations, a value appropriate for the hot Jupiter regime based on analytic modelling by Nichols & Milan (2016).
In the work presented here, simulations were run for an entirely southward IMF orientation, with B sw values ranging from 1-10 6 nT, and with a constant ionospheric Pedersen conductance between 1 -10 4 mho. The assumption of a constant Pedersen conductance provides a reasonable firstorder approximation to the global average. A uniform zero Hall conductance was used throughout, which along with a constant non-zero Pedersen conductance forms the simplest ionospheric conductance model to approximate a realistic magnetospheric configuration. Such a configuration has been used as a standard ionospheric model in many previous MHD simulations (e.g. Fedder & Lyon 1987;Jia et al. 2012a,b).
The simulations were run on a 3D Cartesian computational grid of 256 × 256 × 256 R P . Although BATS-R-US may be run in a time-accurate mode, since this work focuses on hypothetical events, we instead used iterative local time stepping, in which each computation cell takes different time steps and the simulation progresses for a fixed number of iterations to converge on a steady-state solution. At the inner boundary of each simulation the grid was initially highly resolved near the planet with cells of 1/8 R P in size, while the remainder of the grid is incrementally more coarsely resolved moving out further from the planet. This initial resolution is entirely geometric, i.e. it is not based on any physical criteria, but rather is determined based on expectations of where interesting regions of the solution will emerge requiring high resolution. Each simulation was allowed to run for 3,000 iterations before the grid was refined using the adaptive mesh refinement facility within the MHD code. Refinement added 10% more cells in regions of large ∇P and ∇ × B before the simulation resumed running. Refinement based on these criteria was performed every 300 iterations up to a total of 6,000 iterations, giving a final grid containing approximately 15 million cells. After the final refinement the simulation was then allowed to run for a total of 50,000 iterations, by which point the solution had reached an approximate steady-state.
To validate the approach described above, the model was first tested by replicating earlier results of Ridley et al. (2004) and Ridley (2007), the details of which can be found in Appendix A. This work builds on those studies with vastly increased IMF strengths and Pedersen conductances.

Magnetospheric structure
A series of simulation results of the magnetospheric field morphology and plasma density is shown in Fig. 1 for IMF B sw values ranging from 1 -10 5 nT. In each run the Pedersen conductance was initially Σ P = 10 4 mho, representing the highly conductive ionospheres expected at hot Jupiters. The remainder of the input parameters were as stated in the previous section. The magnetic field lines are traced in the Y = 0 plane, with the incident plasma flowing from right to left. Note that in these plots the field line spacing does not necessarily represent magnetic field strength. The lower end of the IMF range represents conditions analogous to the IMF experienced at the Earth. In the region typically associated with hot Jupiters, i.e. 3 -10 stellar radii, an IMF of approximately 10 3 -10 4 nT is expected (Nichols & Milan 2016). Higher IMF strengths have been examined which represent either planets orbiting the star extremely closely (< 3 stellar radii), or planets orbiting stars with exceptionally strong magnetic fields. For instance, at a distance of 10 R * , assuming a predominantly radial field such that B sw ∝ 1/r 2 , a 10 6 nT IMF equates to a star with a surface magnetic field strength of 10 8 nT, approximately three orders of magnitude great than the solar magnetic field. Fig. 1(a) shows that when the planet is exposed to B sw = 1 nT the magnetosphere formed is similar to that at the Earth, i.e. with a clearly visible magnetotail and an apparent magnetosphere on the sub-stellar wind side of the planet. As the IMF is increased B sw = 10 2 nT and B sw = 10 3 nT, Figs. 1(c) and (d) show the lobes of the tail opening as the upstream Alfvén Mach number becomes lower, and the formation of Alfvén wings draped across the planet becomes apparent. At B sw = 10 4 nT, the flow has become sub-Alfvénic, and the Alfvén wings are formed at a large angle from the equatorial plane. As the IMF is increased further to B sw = 10 5 nT the planetary field is dwarfed by the stellar field, and essentially presents no perturbation to the overwhelming IMF. Fig. 2 shows the magnetospheric field line topology and plasma density in close proximity to the planet for the same simulations as Fig. 1. These plots reveal the compression of the sub-stellar wind side magnetosphere due to the IMF pressure as B sw is increased, as well as escape of the dense planetary plasma along open field lines. Note that in Figs. 2(e) and (f) the substellar magnetosphere has collapsed below the 2.5 R p inner boundary of the simulation. Fig. 3 shows a set of ionospheric plots of FAC density and CPCP at a constant Pedersen conductance value of Σ P = 10 4 mho, and with B sw range of 1 − 10 5 nT. In each FAC density plot positive upward current is indicated in red, and negative downward current by blue. The morphology and magnitudes of both quantities within these plots vary as B sw is increased. At low IMF strengths the auroral FAC density is somewhat irregular and diffuse in shape, and is situated at higher colatitudes of ∼ 20−30 • , slightly displaced towards the sub-stellar wind side of the planet. The upward FAC density peaks for B sw = 10 3 nT at 199 µA m −2 . At B sw = 10 4 nT the auroral FAC structure is a narrow oval at ∼ 30 • colatitude, with the peak upward FAC density of 39.97 µA m −2 . As B sw is increased to 10 5 the magnitude of the upward current density falls significantly to 8.67 µA m −2 .

Ionospheric Electrodynamics
A similar saturation and turnover is seen in the CPCP results in Fig. 3(b). From a CPCP potential of ∼ 3 keV for B sw = 1 nT, the potential peaks around 6 keV at B sw = 10 3 nT, before a substantial drop to Φ CPC = 0.457 kV at B sw = 10 5 nT. Saturation of the CPCP, first discussed by Hill et al. (1976), was subsequently tested by Siscoe et al. (2002) against results from MHD simulations, and was modelled in terms of the stellar wind parameters. Various interpretations have been offered for the phenomenon of CPCP saturation. Siscoe et al. (2002) interpreted the saturation as a weakening of the planetary field at the magnetopause due to the field arising from region 1 currents, thus limiting the rate at which reconnection occurs on the sub-stellar wind side of the planet. Alternatively, Kivelson & Ridley (2008) argued that saturation occurs when the solar or stellar wind impedance is greater than the ionospheric impedance, i.e. when the Pedersen conductance Σ P dominates the Alfvén conductance Σ A , causing a partial reflection of Alfvén waves incident on the ionosphere from the solar or stellar wind. The available magnetospheric convection potential is given by (Nichols & Milan 2016) where E sw is the stellar wind motional electric field, and χ is the fraction of the magnetopause standoff distance R mp which constitutes the stellar wind reconnection channel. Observations for Earth determine a value of χ ≈ 0.5 (Milan et al. 2004), and this value is therefore also employed here.
In the reflected Alfvén wave interpretation of Kivelson & Ridley (2008), the electric potential across the ionosphere is given by where the width of the interaction channel, specified by Kivelson & Ridley (2008) as 0.1πR mp is accounted for by the factor ξ = (0.1π χ). Hence, when Σ P Σ A the CPCP tends towards and saturation occurs. For a fixed Pedersen conductance, increased IMF strength leads to a reduced Alfvén conductance, and thus a decrease in CPCP beyond saturation as observed in the results in Fig. 3(b). The saturation effect is also influenced by the decreasing magnetopause standoff distance, but the dominant contributing factor is decreased Alfvén conductance, since the sub-stellar wind side magnetosphere is completely eroded under high IMF strengths. A notable feature of the plots in Figs. 3(a) is the absence of strong region 2 currents. This is an artefact of the MHD   Fig. 1 but for a smaller scale to show the magnetospheric topology in close proximity to the planet code, with several possible causes suggested by Ridley et al. (2001). Since region 2 currents are generated close to the inner boundary of the model, a high resolution is required to produce currents that are even a fraction of the region 1 currents. Increasing the resolution would increase the time taken for simulations to run, and is therefore a trade-off that must be made in consideration of running the models in a timely fashion. Another option which should achieve the same result is to move the inner boundary to a lower value (e.g. from 2.5 R P to 1.5 R P ). As with increasing the resolution, this solution also increases the run time of the simulations. Another possible cause of the weak region 2 currents is that  gradient and curvature drifts of particles at different energies is not addressed by BATS-R-US. The pressure gradient that results from the reconfiguration of the plasma by these drifts may be a source of region 2 currents.

Responses to variable IMF strength and Pedersen conductance
To fully understand the response of the global simulations to the key driving factors, Fig. 4 shows plots of CPCP, total current, maximum ionospheric FAC density, and radio power, versus both B sw and Pedersen conductance. In Figs. 4(a) -(d) the Pedersen conductance was fixed at Σ P = 10 4 mho, since Koskinen et al. (2010) showed that hot Jupiters likely possess highly conductive ionospheres (Σ ≈ 10 4 -10 5 mho), and the modelled parameters are plotted as a function of B sw from 1 -10 6 nT. In Fig. 4(e) -(h) a fixed IMF value of B sw = 10 4 nT was used, and Σ P was varied. The modelled CPCP, shown in Fig. 4(a), initially rises slowly with B sw to a peak of ∼ 6.5 keV at B sw = 10 3 nT before falling away sharply to a value of 4.84 × 10 −2 keV at B sw = 10 6 nT. The analytic expression of equation (16) for CPCP is also plotted for comparison with the results from the SWMF simulations. A discrepancy between the simulation results and the analytic model is apparent: although the general profiles are similar, the SWMF CPCP values are approximately an order of magnitude smaller than the cor-responding analytic values. A point which may be explained by the fact that the analytic model does not account for viscous interactions at the magnetopause boundary (Nichols & Milan 2016). Fig. 4(c) shows the maximum FAC density as a function of B sw . The magnitude of the FAC density in the analytic model is proportional to the CPCP and Pedersen conductance, given by and the full details of this relation can be found in Nichols & Milan (2016). There is reasonable agreement between the SWMF and analytic results, with the saturation and turnover occurring at B sw = 10 3 nT for both, although the simulation results peak at a slighter higher value (∼ 197 µA m −2 ) than the analytic results (∼ 40 µA m −2 ). Integration of FAC density over a hemisphere (using equation 6) demonstrates that the total ionospheric current also slowly increases wih B sw up to a peak of 2.1×10 10 A at B sw = 10 3 nT, with a general profile similar to that of CPCP ( Fig. 4(a)). However, the absolute values of I tot exceed those from the analytic model of Nichols & Milan (2016). Auroral radio power, calculated using equations (7) and (12), is plotted in in Fig.  4(d), and shows a peak value of 4.5×10 14 W at B sw = 10 3 nT, a value approximately four orders of magnitude greater than equivalent peak ECMI emission from Jupiter's aurora.  Cross-polar cap potential, (f) total current, (g) maximum fieldaligned current density, and (h) auroral radio power as functions of Pedersen conductance for an exoplanet exposed to B sw = 10 4 nT. Blue lines represent the SWMF results, with diamonds denoting the values are which MHD simulations were conducted. Black lines show the analytic results using the model of Nichols & Milan (2016).

Figs. 4(e) -(h)
show the same parameters plotted as a function of Pedersen conductance for a fixed IMF strength of B sw = 10 4 nT. The results show that the total maximum FAC density, total current, and radio power are all virtually independent of Pedersen conductance for the SWMF simulations, since the low Alfvén conductance implied by the high IMF strength means that the condition for saturation described in Section 4.2, namely when Σ P Σ A , is now satisfied for low Σ P values.

Variable orbital distance
Planets at different orbital distances from the host star are not only subject to varying IMF strengths, but also varying stellar wind velocity, density, and Pedersen conductance. In this section the effects on the CPCP, FACs, and radio power are investigated as functions of orbital distance. For a Sun- 10 0 10 1 10 2 10 3 10 4 10 5 10 6 Bsw /nT like star, Nichols & Milan (2016) calculated analytically how stellar wind parameters vary with orbital distance. Using that work as a reference, spot values of stellar wind velocity, density, IMF strength, and Pedersen conductance were taken at five orbital distances from 2 -126 R * , and used as inputs for five SWMF simulations. The input parameter values used for each run in the simulation set are summarised in Table 1, and Fig. 5 shows the Pedersen conductance and IMF values used, in relation to the cuts in the parameter space representative of the simulation sets in Section 4.3. Fig. 6(a) shows that the CPCP rises from ∼ 3 kV in the nominal hot Jupiter orbital distance region, to ∼ 106 kV at an orbital distance of 126 R * . Total FAC (Fig. 6(b)) falls with increasing orbital distance from a value of ∼ 5 × 10 10 A at d = 2 R * to ∼ 2 × 10 9 A at d = 126 R * , and a similar trend is observed in the results for maximum FAC density (Fig.  6(c)). Figs. 6(a) -(c) show a good agreement between the SWMF results and the analytic results of Nichols & Milan (2016). In Fig. 6(d) auroral radio power is shown along with analytic results of Nichols & Milan (2016), and the Radiometric Bode's law. The dashed line shows a least-squares polynomial fit to the SWMF results, which yields the relation P r ∝ d −1.398 . Nichols & Milan (2016), using a Parker spiral IMF, found that the power varies as P r ∝ d −5/2 in the inner orbtial distance range, i.e. before the notch in Fig.  6(d), and as P r ∝ d −5/4 in the outer orbital distance range, where the resultant IMF is dominated by the perpendicular field component. The relation found in this study of P r ∝ d −1.398 therefore lies between the two power laws determined by Nichols & Milan (2016). The radio powers in the hot Jupiter orbital range (3 -10 R * ) of ∼ 10 15 W are commensurate with the peak radio powers seen in the results of Fig. 4. Finally, Fig. 6(e) plots the spectral flux density F r calculated using equation (13) at a distance of s = 15 pc, a value chosen since it is apparent that emission from planets significantly beyond this distance would be below the detection threshold of currently available radio telescopes. Horizontal lines in Fig. 6(e) indicate the sensitivities of the Murchison Widefield Array (MWA), the Low-Frequency Array (LOFAR), the Very Large Array (VLA), and the Square Kilometre Array (SKA). The results show that radio flux from hot Jupiters located within 15 pc of the Solar system should be detectable with both the VLA, and in the future with the SKA.

DISCUSSION AND CONCLUSIONS
The model used in the study contains some inherent limitations due to simplifications made to the realistic dynamics of star-hot Jupiter interactions. The upper atmospheres of hot Jupiters undergo intensive escape due to ionization and radiation heating driven by stellar X-ray and EUV radiation, giving rise to the so-called planetary wind (e.g. Yelle 2004;Muñoz 2007;Murray-Clay et al. 2009;Koskinen et al. 2013). This planetary wind is additionally shaped by the gravitational interaction between the star and planet. As can be seen from equations (1) -(5), gravitational effects and ionisation due to heating and radiation are not included in the MHD model used in this study. Atmospheric escape in the form of a planetary wind would add an additional pressure outward from the planet which would expand the magnetosphere, and we expect this to increase the strength of the radio emission due to an increase in the size of the stellar wind reconnection channel. This study is a first attempt to establish the M-I coupling dynamics using a 3D MHD model in the framework of existing analytic studies, and as such is primarily intended to investigate the broadbrush effect of IMF strength and Pedersen conductance on the magnetospheric FACs and auroral radio emission. Attempting to incorporate the phenomena of gravitational and ionisation effects immediately may introduce additional unconstrained free parameters which could obscure the intention of this study, namely to isolate the response of the model outputs to IMF strength and ionospheric Pedersen conductance. A more complete future study should develop this model to incorporate the gravitational and ionization effects on the star-planet interaction and resulting FACs and radio emission.
One of the most notable features of the results shown in Fig. 4 is the approximately order of magnitude discrepancy between the numerical SWMF results and analytic model for the CPCP. A possible explanation for this is suggested by Ridley et al. (2004), who remark that the magnetospheres studied in MHD simulations are, by nature, MHD magnetospheres. Despite accurately depicting the general shape of the pressure distribution, the magnitude is typically underestimated by approximately an order of magnitude at the inner magnetosphere. This underestimation is a result of the lack of energy discretization in the MHD simulation, meaning that the code is unable to model high energy particles. As the FAC densities in the simulations presented here are substantially higher than those typically encountered in similar SWMF modelling of Earth magnetosphere, the associated electrons also have correspondingly higher energies. Hence, an inability to model high energy electron may have more impact to this study than to studies involving less energetic particles at Earth, and thus provide a plausible explanation for the order of magnitude discrepancy seen between numerical and analytic results.
The saturation of field-aligned currents and radio power at relatively low ionospheric Pedersen conductances implies that variations in this parameter may be largely inconsequential at hot Jupiters, where anticipated conductance values are of the order 10 4 -10 5 mho. Analytic modelling finds that Pedersen conductance affects the point at which CPCP, and therefore radio power, saturates, but this effect is not apparent in the results from SWMF simulations, where ra-dio power peaks at a point which appears to be independent of Pedersen conductance. Fig. 4 shows that as Pedersen conductance is increased, CPCP falls while FAC density remains largely constant for the SWMF results. This reinforces previous findings that the magnetosphere acts as neither a voltage nor current generator (Fedder & Lyon 1987;Ridley et al. 2004). By considering Ohm's law, J = σE, if the magnetosphere is acting as a current generator, where J is constant, then increasing σ would result in a linear increase of E. However, Fig. 4 shows a constant current but non-linear decreasing CPCP. One the other hand, if the magnetosphere acts as a voltage generator then increasing σ would correspond to a linear increase of J for a constant E. Neither of these trends are seen in the SWMF results, and in the analytic results both current and potential are found to change simultaneously with Σ P . The results shown in Fig. 4 represent spot values at horizontal and vertical cuts along a B sw -Σ P plane, and exhibit differences between the analytic and numerical MHD values. However, the values plotted in Fig. 6 for a realistic variation with orbital distance actually agree remarkably well.
Various scaling laws exist which approximate the planetary magnetic dynamo performance to provide estimations of the magnetic field strengths at hot Jupiters. For example Griessmeier et al. (2005Griessmeier et al. ( , 2007 predict magnetic moments of hot Jupiters to be approximately 10% of the Jovian field strength. However, a recent study by Cauley et al. (2019) find evidence for hot Jupiter surface magnetic field strengths approximately 10 -100 times larger than those predicted by scaling laws. Therefore, while magnetic field strengths at hot Jupiters may be a fraction of the Jovian value, it is also possible that they exceed the Jovian value. In light of the present high degree of uncertainty in the field regarding hot Jupiter magnetic field strengths determined either from scaling laws or from observations of star-planet interaction, the jovian value is employed in this study as a reasonable benchmark alongside which the jovian values of other presently unconstrained hot Jupiter parameters (i.e. source plasma population number density and thermal energy) may also be employed for consistency.
Note that we have not considered here hydrodynamic outflow owing to the strong irradiation from the host star or stellar-planetary tidal interaction, and as such our work should be considered to be a first step toward a selfconsistent picture of the magnetospheric dynamics of a hot Jupiter. The escaping planetary wind at hot Jupiters should deform the field lines from their dipolar topology, and under certain conditions may distort the magnetosphere into a magnetodisc configuration (Khodachenko et al. 2015). At Jupiter the modification of the planetary magnetosphere into a disc-like topology affects the FAC system and the location and size of the main auroral oval (Nichols 2011;Nichols et al. 2015). Similar effects may occur at hot Jupiters, although the presence of a magnetodisc may not be typical for each hot Jupiter, and the effect depends on a variety of parameters which are beyond the scope of this present paper and would require a dedicated future study.
The use of a constant Pedersen conductance in this work serves as a first approximation to a realistic ionosphere, but future work should examine more plausible conductance models. Exoplanets in close orbit around the host star would likely be tidally locked, meaning that one side of the planet permanently faces the star, and is subject to intense ionising stellar radiation. This could result in an asymmetric Pedersen conductance pattern, with an ionosphere which is highly conductive on the sub-stellar wind side of the planet and has a low conductance on the opposite side. Particle precipitation from the auroral currents may also further ionise the atmosphere, amplifying the conductance. Such a self-consistent ionospheric model may modify the result presented in this paper, but the general findings would not be expected to alter significantly. This work may also be extended to incorporate planetary rotation. Since the simulations presented here were run in a time-independent mode, planetary rotation was not a factor, but with simple modification simulations could be run in a time-accurate mode to investigate the effects of rotation on the field-aligned currents and radio power. In time-accurate mode it would also be possible to examine different stellar wind configurations from the entirely southward B z IMF in this work. Hence a more realistic Parker spiral type IMF, or a north-south switching IMF could both be implemented. A future study to explore planetary magnetic field strengths greater the jovian value employed here is also warranted. The trade-off to consider in doing so is the increase in run-time for the SWMF.
Ultimately, this work is intended to guide observations regarding the feasibility of detecting auroral radio emission from exoplanets. The maximum predicted radio powers in this study of ∼ 10 14 -10 15 W are consistent with the findings of Nichols & Milan (2016), but, in the hot Jupiter regime, are lower than predicted by studies employing the RBL (e.g. Zarka et al. 2001), potentially explaining the lack of detection to date. A particular finding in this work is that such intense radio emission may only occur for a planet exposed to a narrow range of stellar magnetic field strengths. However, this assumes other stellar wind parameters are unchanging, and the simulations for a range of orbital distances (Fig. 6), taking into account the variation of these other stellar wind parameters, show a relation between radio power and orbital distance similar to that found in the analytic work of Nichols & Milan (2016). The flux density calculations in this study suggest that auroral emission directly from hot Jupiters in the local stellar neighbourhood (≤ 15 pc) may be presently detectable with telescopes such as the VLA, and in the near future with the SKA. Separating an exoplanetary radio signal from background noise may prove challenging, although modulation of the signal at planetary orbital periods could allow light curves to be folded at that period to aid identification of radio signals. Ridley et al. (2004) and Ridley (2007) investigated the influence of variable IMF strength and Pedersen conductance on ionospheric cross-polar cap potential and field-aligned currents for Earth. Since this paper investigates similar effects, but at much larger scales and magnitudes, in order to test the method used and to validate the results of this work, the key results from Ridley et al. (2004) and Ridley (2007) were first reproduced and compared with results from the original studies.

Studies by
Simulations using the SWMF were initialised with the terrestrial radius, dipole field strength, and plasma density, as per Ridley et al. (2004); Ridley (2007). Adaptive mesh refinement of the computational grid was performed in a manner similar to that described by Ridley (2007). Fig. A1 shows the results for the ionospheric field-aligned current density and cross-polar cap potential as functions of Pedersen conductance, which compare closely with the results of Ridley et al. (2004) (Fig. A2). Fig. A3 plots cross-polar cap potential as a function of IMF strength, and again the results are in reasonable agreement with those of Ridley (2007) (Fig.  A4), thus validating the technique employed in this study. This paper has been typeset from a T E X/L A T E X file prepared by the author.