# The space weather around the exoplanet GJ 436b. II. Stellar wind-exoplanet interactions

A. A. Vidotto, V. Bourrier, R. Fares, S. Bellotti, J. F. Donati, P. Petit, G. A.

J. Hussain, Julien Morin

► To cite this version:

A. A. Vidotto, V. Bourrier, R. Fares, S. Bellotti, J. F. Donati, et al. The space weather around the exoplanet GJ 436b. II. Stellar wind-exoplanet interactions. Astronomy and Astrophysics - A&A, 2023, 678, 10.1051/0004-6361/202347237. insu-04473152

HAL Id: insu-04473152

https://insu.hal.science/insu-04473152

Submitted on 23 Feb 2024

**HAL** is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire **HAL**, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# The space weather around the exoplanet GJ 436b

## II. Stellar wind–exoplanet interactions

A. A. Vidotto<fnref n="1" />, V. Bourrier<fnref n="2" />, R. Fares<fnref n="3" />, S. Bellotti<fnref n="1" />,<fnref n="4" />,<fnref n="5" />, J. F. Donati<fnref n="5" />, P. Petit<fnref n="5" />, G. A. J. Hussain<fnref n="4" />, and J. Morin<fnref n="6" />

Received 20 June 2023 / Accepted 23 August 2023

### ABSTRACT

The M dwarf star GJ 436 hosts a warm-Neptune that is losing a substantial amount of atmosphere, which is shaped by the interactions with the wind of the host star. The stellar wind is formed by particles and magnetic fields that shape the exo-space weather around the exoplanet GJ 436b. Here, we use the recently published magnetic map of GJ 436 to model its 3D Alfvén-wave-driven wind. We compared our results with previous transmission spectroscopic models and measurements of non-thermal velocities at the transition region of GJ 436; our models indicate that the wind of GJ 436 is powered by a smaller flux of Alfvén waves than that powering the wind of the Sun. This suggests that the canonical flux of Alfvén waves assumed in solar wind models might not be applicable to the winds of old M dwarf stars. Compared to the solar wind, GJ 436's wind has a weaker acceleration and an extended sub-Alfvénic region. This is important because it places the orbit of GJ 436b inside the region dominated by the stellar magnetic field (i.e. inside the Alfvén surface). Due to the sub-Alfvénic motion of the planet through the stellar wind, magnetohydrodynamic waves and particles released in reconnection events can travel along the magnetic field lines towards the star, which could power the anomalous ultraviolet flare distribution recently observed in the system. For an assumed planetary magnetic field of $B_p \approx 2$ G, we derive the power released by stellar wind–planet interactions as $\mathcal{P} \sim 10^{22}-10^{23}$ erg s$^{-1}$, which is consistent with the upper limit of $10^{26}$ erg s$^{-1}$ derived from ultraviolet lines. We further highlight that, because star–planet interactions depend on stellar wind properties, observations that probe these interactions and the magnetic map used in 3D stellar wind simulations should be contemporaneous for deriving realistic results.

**Key words.** stars: winds, outflows – stars: individual: GJ436 – planetary systems – planet-star interactions

## 1. Introduction

The evolution of planets is shaped by billions of years of interactions with their host stars and with neighbouring planets. In addition to the central gravitational force of the host star, there are several types of interactions taking place between a star and an exoplanet (Vidotto 2020), some of which are not present, or simply very weak, in the present-day Solar System planets.

In the case of the GJ 436 system, the focus of the current study, the warm Neptune GJ 436b experiences a strong interaction with its host star due to the combined effects of high-energy stellar irradiation that is deposited in the planetary atmosphere and a subsequent interaction of this atmosphere with the stellar wind (Bourrier et al. 2016). Spectroscopic transit observations in Lyman-α (Kulow et al. 2014; Ehrenreich et al. 2015; Lavie et al. 2017; Dos Santos et al. 2019) revealed that GJ 436b is enshrouded by a giant cloud of escaping atmosphere that trails GJ 436b along its orbit. Three-dimensional models suggest that this trailing structure takes the form of a giant comet-like tail (Bourrier et al. 2016; Shaikhislamov et al. 2018; Khodachenko et al. 2019; Villarreal D'Angelo et al. 2021). By

modelling the observations of spectroscopic transits of GJ 436b with atmospheric escape models that account for atmospheric interaction with the stellar wind, local properties of the stellar wind at the orbital distance of the planet, such as the wind speed and density, can be estimated; they can then be used to infer global properties of stellar winds, such as their mass-loss rates (Vidotto & Bourrier 2017). All these quantities shape the space weather around a planet.

Indeed, one of the major strengths in observing and modelling star–planet interactions in general is that they can provide alternative ways to physically characterise planetary systems. Stellar wind mass-loss rates, which are difficult to constrain in the case of cool dwarf stars (Vidotto 2021), can also be estimated by investigating star–planet interactions through planet-induced radio emission (Kavanagh et al. 2021; Pineda & Villadsen 2023) or auroral radio emission from the exoplanet (Vidotto et al. 2019). Investigating star–planet interactions through chromospheric hot spots associated with anomalous stellar activity (e.g. Shkolnik et al. 2005; Cauley et al. 2019) or exoplanetary radio emission (e.g. Farrell et al. 1999; Turner et al. 2021) allows one to probe planetary magnetism, a quantity that is believed to be

<footnote n="1">¹ Leiden Observatory, Leiden University, PO Box 9513, 2300 RA, Leiden, The Netherlands
e-mail: vidotto@strw.leidenuniv.nl</footnote>

<footnote n="2">² Observatoire Astronomique de l'Université de Genève, Chemin Pegasi 51b, 1290 Versoix, Switzerland</footnote>

<footnote n="3">³ Department of Physics, College of Science, United Arab Emirates University, PO Box No. 15551, Al Ain, UAE
e-mail: rim fares@uaeu.ac.ae</footnote>

<footnote n="4">⁴ Science Division, Directorate of Science, European Space Research and Technology Centre (ESA/ESTEC), Keplerlaan 1, 2201 AZ, Noordwijk, The Netherlands</footnote>

<footnote n="5">⁵ Institut de Recherche en Astrophysique et Planétologie, Université de Toulouse, CNRS, IRAP/UMR 5277, 14 avenue Édouard Belin, 31400 Toulouse, France</footnote>

<footnote n="6">⁶ Laboratoire Univers et Particules de Montpellier, Université de Montpellier, CNRS, 34095 Montpellier, France</footnote>
# The space weather around the exoplanet GJ 436b

## II. Stellar wind–exoplanet interactions

A. A. Vidotto<fnref n="1" />, V. Bourrier<fnref n="2" />, R. Fares<fnref n="3" />, S. Bellotti<fnref n="1" />,<fnref n="4" />,<fnref n="5" />, J. F. Donati<fnref n="5" />, P. Petit<fnref n="5" />, G. A. J. Hussain<fnref n="4" />, and J. Morin<fnref n="6" />

Received 20 June 2023 / Accepted 23 August 2023

## ABSTRACT

The M dwarf star GJ 436 hosts a warm-Neptune that is losing a substantial amount of atmosphere, which is shaped by the interactions with the wind of the host star. The stellar wind is formed by particles and magnetic fields that shape the exo-space weather around the exoplanet GJ 436b. Here, we use the recently published magnetic map of GJ 436 to model its 3D Alfvén-wave-driven wind. We compared our results with previous transmission spectroscopic models and measurements of non-thermal velocities at the transition region of GJ 436; our models indicate that the wind of GJ 436 is powered by a smaller flux of Alfvén waves than that powering the wind of the Sun. This suggests that the canonical flux of Alfvén waves assumed in solar wind models might not be applicable to the winds of old M dwarf stars. Compared to the solar wind, GJ 436's wind has a weaker acceleration and an extended sub-Alfvénic region. This is important because it places the orbit of GJ 436b inside the region dominated by the stellar magnetic field (i.e. inside the Alfvén surface). Due to the sub-Alfvénic motion of the planet through the stellar wind, magnetohydrodynamic waves and particles released in reconnection events can travel along the magnetic field lines towards the star, which could power the anomalous ultraviolet flare distribution recently observed in the system. For an assumed planetary magnetic field of $B_p \approx 2$ G, we derive the power released by stellar wind–planet interactions as $\mathcal{P} \sim 10^{22}-10^{23}$ erg s$^{-1}$, which is consistent with the upper limit of $10^{26}$ erg s$^{-1}$ derived from ultraviolet lines. We further highlight that, because star–planet interactions depend on stellar wind properties, observations that probe these interactions and the magnetic map used in 3D stellar wind simulations should be contemporaneous for deriving realistic results.

**Key words.** stars: winds, outflows – stars: individual: GJ436 – planetary systems – planet-star interactions

## 1. Introduction

The evolution of planets is shaped by billions of years of interactions with their host stars and with neighbouring planets. In addition to the central gravitational force of the host star, there are several types of interactions taking place between a star and an exoplanet (Vidotto 2020), some of which are not present, or simply very weak, in the present-day Solar System planets.

In the case of the GJ 436 system, the focus of the current study, the warm Neptune GJ 436b experiences a strong interaction with its host star due to the combined effects of high-energy stellar irradiation that is deposited in the planetary atmosphere and a subsequent interaction of this atmosphere with the stellar wind (Bourrier et al. 2016). Spectroscopic transit observations in Lyman-α (Kulow et al. 2014; Ehrenreich et al. 2015; Lavie et al. 2017; Dos Santos et al. 2019) revealed that GJ 436b is enshrouded by a giant cloud of escaping atmosphere that trails GJ 436b along its orbit. Three-dimensional models suggest that this trailing structure takes the form of a giant comet-like tail (Bourrier et al. 2016; Shaikhislamov et al. 2018; Khodachenko et al. 2019; Villarreal D'Angelo et al. 2021). By

modelling the observations of spectroscopic transits of GJ 436b with atmospheric escape models that account for atmospheric interaction with the stellar wind, local properties of the stellar wind at the orbital distance of the planet, such as the wind speed and density, can be estimated; they can then be used to infer global properties of stellar winds, such as their mass-loss rates (Vidotto & Bourrier 2017). All these quantities shape the space weather around a planet.

Indeed, one of the major strengths in observing and modelling star–planet interactions in general is that they can provide alternative ways to physically characterise planetary systems. Stellar wind mass-loss rates, which are difficult to constrain in the case of cool dwarf stars (Vidotto 2021), can also be estimated by investigating star–planet interactions through planet-induced radio emission (Kavanagh et al. 2021; Pineda & Villadsen 2023) or auroral radio emission from the exoplanet (Vidotto et al. 2019). Investigating star–planet interactions through chromospheric hot spots associated with anomalous stellar activity (e.g. Shkolnik et al. 2005; Cauley et al. 2019) or exoplanetary radio emission (e.g. Farrell et al. 1999; Turner et al. 2021) allows one to probe planetary magnetism, a quantity that is believed to be

<footnote n="1">¹ Leiden Observatory, Leiden University, PO Box 9513, 2300 RA, Leiden, The Netherlands
e-mail: vidotto@strw.leidenuniv.nl</footnote>

<footnote n="2">² Observatoire Astronomique de l'Université de Genève, Chemin Pegasi 51b, 1290 Versoix, Switzerland</footnote>

<footnote n="3">³ Department of Physics, College of Science, United Arab Emirates University, PO Box No. 15551, Al Ain, UAE
e-mail: rim fares@uaeu.ac.ae</footnote>

<footnote n="4">⁴ Science Division, Directorate of Science, European Space Research and Technology Centre (ESA/ESTEC), Keplerlaan 1, 2201 AZ, Noordwijk, The Netherlands</footnote>

<footnote n="5">⁵ Institut de Recherche en Astrophysique et Planétologie, Université de Toulouse, CNRS, IRAP/UMR 5277, 14 avenue Édouard Belin, 31400 Toulouse, France</footnote>

<footnote n="6">⁶ Laboratoire Univers et Particules de Montpellier, Université de Montpellier, CNRS, 34095 Montpellier, France</footnote>
very important in the context of planetary habitability but is still
rather poorly known in planets outside the Solar System.

In this work we focus on star–planet interactions that are
mediated by the stellar wind and its embedded magnetic field. In
the literature, this type of interaction is referred to as ‘star–planet
(electro)magnetic interactions’, ‘stellar wind–planet interactions’
or ‘electrodynamic star–planet interactions’, but they are only
different names of the same process. In the context of star–planet
interactions mediated by the magnetised host star’s wind, differ-
ent physical processes take place depending on the wind regime
in which the planet orbits.

In the regime of sub-Alfvénic motion of a planet through
the wind of its host star, a magnetic ‘connectivity’ between the
planet and the star can take place. In this regime, the orbiting
planet can trigger magnetohydrodynamic (MHD) waves that can
travel towards the star, or alternatively (and even possibly con-
currently) it can also trigger magnetic reconnection events in the
case of a magnetised planet. Both of these mechanisms release
energy that is radiated away and could be detected in the form of,
for example, chromospheric hot spots (Lanza 2009; Strugarek
et al. 2022), planet-induced coronal radio emission (Saur
et al. 2013; Pineda & Villadsen 2023), anomalous flare events
(Fischer & Saur 2019), and even potentially the lack of energetic
flare events, as recently proposed by Loyd et al. (2023) in the
context of GJ 436. We note that all these signatures originate on
or near the host star, and not on the planet.

The inverse takes place in the regime of super-Alfvénic
motion of a planet through the wind of its host star. In this
regime, the interaction between the planet and the stellar wind
can lead to radio emission originating from the planet itself
(more precisely, from its magnetosphere). This is because, in the
super-Alfvénic regime, the magnetic connectivity between the
star and planet can no longer form. This is the situation currently
experienced by the Solar System planets, where we see that the
planetary auroral radio emission is correlated with the solar wind
energy that is dissipated in the magnetosphere of Solar System
planets (e.g. Farrell et al. 1999). Despite many attempts to detect
such signatures in exoplanets, it was only recently that a poten-
tial detection of radio emission from a hot Jupiter was reported
(Turner et al. 2021, see also Elekes & Saur 2023 for models of
these observations).

Fundamental to the nature of stellar wind–planet interactions
is thus the regime of the orbital motion of the planet: sub- or
super-Alfvénic. To determine the regime, we need to obtain more
realistic stellar wind properties, such as densities, velocities, and
magnetism. These stellar wind quantities are also what shape
the space weather of planets and are thus also relevant for, for
example, understanding the flux of energetic particles impacting
planetary atmospheres (e.g. Mesquita et al. 2021; Rodgers-Lee
et al. 2023). In the case of GJ 436, its wind has been investi-
gated in several works (Vidotto & Bourrier 2017; Mesquita &
Vidotto 2020; Mesquita et al. 2021), but it is only recently that
the surface magnetic field of the host star GJ 436 was mapped
(Bellotti et al. 2023); this allows us to now model, for the first
time, the 3D structure of the stellar wind, including the observed
surface magnetism. Magnetic fields play an important role in the
winds of cool dwarf stars, as they are responsible for the mech-
anisms that heat and accelerate these winds. Section 2 presents
our 3D MHD simulations of the wind of GJ 436. We then use the
results of our wind models to compute the space weather con-
ditions (Sect. 3) along the highly misaligned orbit of GJ 436b.
GJ 436b is a Neptune-mass planet orbiting very close to its
slowly rotating M dwarf host. It has a puzzling observed orbital
architecture – a non-circular orbit with an eccentricity of 0.152

<figure><img src="image_1.png" /><figcaption>**Fig. 1.** Observationally reconstructed radial magnetic field of GJ 436 from Bellotti et al. (2023) using Zeeman-Doppler imaging. This map is incorporated at the boundary of our stellar wind simulations.</figcaption></figure>

(Trifonov et al. 2018) and strong orbital misalignment (Bourrier et al. 2018, 2022) – that is likely the result of a delayed migration following Kozai-Lidov interactions with an outer massive companion (Bourrier et al. 2018; Attia et al. 2021). In Sect. 4 we then compute the strengths of star–planet interactions. In the final section of this paper (Sect. 5), we discuss how our computations compare to the recent observations of Loyd et al. (2023) and present our concluding thoughts.

## 2. The wind of GJ 436

The environment surrounding GJ 436b consists of particles and
magnetic fields that form the wind of the host star. Here, we char-
acterise this environment by performing data-driven 3D MHD
simulations of the stellar wind. Because winds of cool dwarf
stars are ultimately driven by stellar magnetism, at the bound-
ary of our model (which is located close to the stellar surface, at
the transition region), we include the reconstructed surface field
shown in Fig. 1 (Bellotti et al. 2023). This map was reconstructed
using the Zeeman-Doppler imaging technique (Donati & Brown
1997; Donati & Landstreet 2009) and, thus, represents the large-
scale field of GJ 436 in March-June 2016. For this simulation, we
used the Alfvén Wave Solar Model (AWSOM; van der Holst et al.
2014) implemented in the numerical code BATS-R-US (Powell
et al. 1999; Toth et al. 2012). We refer the reader to van der Holst
et al. (2014) for further details of the model. Below we describe
the main input parameters adopted in the simulations of the wind
of GJ 436.

In addition to the stellar surface magnetic map, other model
inputs are the stellar properties (mass, radius and rotation period;
cf. Table 1), the density and temperature at the transition region,
and the Alfvén wave properties (boundary wave flux and the
correlation length of the Alfvénic wave turbulence). Of these
parameters, the least constrained ones are the wave properties,
as we discuss below.

The inner boundary of our model is at the transition region,
and values considered for solar wind models (van der Holst
et al. 2014; Oran et al. 2017) for ion/electron density ($n_{tr}$ =
$5 \times 10^{11}$ cm$^{-3}$), temperature ($T_{tr}$ = $5 \times 10^{4}$ K), and, conse-
quently, gas density ($P_{tr}$ ≈ 1.4 dyn cm$^{-2}$), were adopted. Our
assumed pressure is not too dissimilar to the value assumed
in the chromospheric models of the M dwarf star GJ832 by
Fontenla et al. (2016), who used $P_{tr}$ = 0.8 dyn cm$^{-2}$. For a tem-
perature of $5 \times 10^{4}$ K, this pressure value translates to a density
of $10^{11}$ cm$^{-3}$, which is also similar to density values found at this
temperature from the chromosphere+wind models of Sakaue &
Shibata (2021). However, there are also models in the literature
that adopt lower pressure values of $P_{tr}$ = 0.02–0.04 dyn cm$^{-2}$
very important in the context of planetary habitability but is still
rather poorly known in planets outside the Solar System.

In this work we focus on star–planet interactions that are
mediated by the stellar wind and its embedded magnetic field. In
the literature, this type of interaction is referred to as ‘star–planet
(electro)magnetic interactions’, ‘stellar wind–planet interactions’
or ‘electrodynamic star–planet interactions’, but they are only
different names of the same process. In the context of star–planet
interactions mediated by the magnetised host star’s wind, differ-
ent physical processes take place depending on the wind regime
in which the planet orbits.

In the regime of sub-Alfvénic motion of a planet through
the wind of its host star, a magnetic ‘connectivity’ between the
planet and the star can take place. In this regime, the orbiting
planet can trigger magnetohydrodynamic (MHD) waves that can
travel towards the star, or alternatively (and even possibly con-
currently) it can also trigger magnetic reconnection events in the
case of a magnetised planet. Both of these mechanisms release
energy that is radiated away and could be detected in the form of,
for example, chromospheric hot spots (Lanza 2009; Strugarek
et al. 2022), planet-induced coronal radio emission (Saur
et al. 2013; Pineda & Villadsen 2023), anomalous flare events
(Fischer & Saur 2019), and even potentially the lack of energetic
flare events, as recently proposed by Loyd et al. (2023) in the
context of GJ 436. We note that all these signatures originate on
or near the host star, and not on the planet.

The inverse takes place in the regime of super-Alfvénic
motion of a planet through the wind of its host star. In this
regime, the interaction between the planet and the stellar wind
can lead to radio emission originating from the planet itself
(more precisely, from its magnetosphere). This is because, in the
super-Alfvénic regime, the magnetic connectivity between the
star and planet can no longer form. This is the situation currently
experienced by the Solar System planets, where we see that the
planetary auroral radio emission is correlated with the solar wind
energy that is dissipated in the magnetosphere of Solar System
planets (e.g. Farrell et al. 1999). Despite many attempts to detect
such signatures in exoplanets, it was only recently that a poten-
tial detection of radio emission from a hot Jupiter was reported
(Turner et al. 2021, see also Elekes & Saur 2023 for models of
these observations).

Fundamental to the nature of stellar wind–planet interactions
is thus the regime of the orbital motion of the planet: sub- or
super-Alfvénic. To determine the regime, we need to obtain more
realistic stellar wind properties, such as densities, velocities, and
magnetism. These stellar wind quantities are also what shape
the space weather of planets and are thus also relevant for, for
example, understanding the flux of energetic particles impacting
planetary atmospheres (e.g. Mesquita et al. 2021; Rodgers-Lee
et al. 2023). In the case of GJ 436, its wind has been investi-
gated in several works (Vidotto & Bourrier 2017; Mesquita &
Vidotto 2020; Mesquita et al. 2021), but it is only recently that
the surface magnetic field of the host star GJ 436 was mapped
(Bellotti et al. 2023); this allows us to now model, for the first
time, the 3D structure of the stellar wind, including the observed
surface magnetism. Magnetic fields play an important role in the
winds of cool dwarf stars, as they are responsible for the mech-
anisms that heat and accelerate these winds. Section 2 presents
our 3D MHD simulations of the wind of GJ 436. We then use the
results of our wind models to compute the space weather con-
ditions (Sect. 3) along the highly misaligned orbit of GJ 436b.
GJ 436b is a Neptune-mass planet orbiting very close to its
slowly rotating M dwarf host. It has a puzzling observed orbital
architecture – a non-circular orbit with an eccentricity of 0.152

<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Observationally reconstructed radial magnetic field of GJ 436 from Bellotti et al. (2023) using Zeeman-Doppler imaging. This map is incorporated at the boundary of our stellar wind simulations.</figcaption></figure>

(Trifonov et al. 2018) and strong orbital misalignment (Bourrier et al. 2018, 2022) – that is likely the result of a delayed migration following Kozai-Lidov interactions with an outer massive companion (Bourrier et al. 2018; Attia et al. 2021). In Sect. 4 we then compute the strengths of star–planet interactions. In the final section of this paper (Sect. 5), we discuss how our computations compare to the recent observations of Loyd et al. (2023) and present our concluding thoughts.

## 2. The wind of GJ 436

The environment surrounding GJ 436b consists of particles and
magnetic fields that form the wind of the host star. Here, we char-
acterise this environment by performing data-driven 3D MHD
simulations of the stellar wind. Because winds of cool dwarf
stars are ultimately driven by stellar magnetism, at the bound-
ary of our model (which is located close to the stellar surface, at
the transition region), we include the reconstructed surface field
shown in Fig. 1 (Bellotti et al. 2023). This map was reconstructed
using the Zeeman-Doppler imaging technique (Donati & Brown
1997; Donati & Landstreet 2009) and, thus, represents the large-
scale field of GJ 436 in March-June 2016. For this simulation, we
used the Alfvén Wave Solar Model (AWSOM; van der Holst et al.
2014) implemented in the numerical code BATS-R-US (Powell
et al. 1999; Toth et al. 2012). We refer the reader to van der Holst
et al. (2014) for further details of the model. Below we describe
the main input parameters adopted in the simulations of the wind
of GJ 436.

In addition to the stellar surface magnetic map, other model
inputs are the stellar properties (mass, radius and rotation period;
cf. Table 1), the density and temperature at the transition region,
and the Alfvén wave properties (boundary wave flux and the
correlation length of the Alfvénic wave turbulence). Of these
parameters, the least constrained ones are the wave properties,
as we discuss below.

The inner boundary of our model is at the transition region,
and values considered for solar wind models (van der Holst
et al. 2014; Oran et al. 2017) for ion/electron density ($n_{tr}$ =
$5 \times 10^{11}$ cm$^{-3}$), temperature ($T_{tr}$ = $5 \times 10^{4}$ K), and, conse-
quently, gas density ($P_{tr}$ ≈ 1.4 dyn cm$^{-2}$), were adopted. Our
assumed pressure is not too dissimilar to the value assumed
in the chromospheric models of the M dwarf star GJ832 by
Fontenla et al. (2016), who used $P_{tr}$ = 0.8 dyn cm$^{-2}$. For a tem-
perature of $5 \times 10^{4}$ K, this pressure value translates to a density
of $10^{11}$ cm$^{-3}$, which is also similar to density values found at this
temperature from the chromosphere+wind models of Sakaue &
Shibata (2021). However, there are also models in the literature
that adopt lower pressure values of $P_{tr}$ = 0.02–0.04 dyn cm$^{-2}$
<table><caption>**Table 1.** Parameters adopted in our models for GJ 436 and GJ 436b.</caption><thead><tr><th>Parameter (Star)</th><th>Value</th><th>Units</th><th>Reference</th></tr></thead><tbody><tr><td>Mass (M<sub>★</sub>)</td><td>≈0.45</td><td>M<sub>☉</sub></td><td>Maxted et al. (2022)</td></tr><tr><td>Radius (R<sub>★</sub>)</td><td>≈0.42</td><td>R<sub>☉</sub></td><td>Maxted et al. (2022)</td></tr><tr><td>Rotation period (P<sub>rot</sub>)</td><td>≈44.1</td><td>days</td><td>Bourrier et al. (2018)</td></tr><tr><td>Stellar magnetism</td><td>cf. Fig. 1</td><td></td><td>Bellotti et al. (2023)</td></tr><tr><td>Age</td><td>6<sup>+4</sup><sub>-5</sub></td><td>Gyr</td><td>Torres (2007)</td></tr></tbody></table><table><thead><tr><th>Parameter (Planet)</th><th>Value</th><th>Units</th><th>Reference</th></tr></thead><tbody><tr><td>Mass (M<sub>pl</sub>)</td><td>21.72</td><td>M<sub>⊕</sub></td><td>Maxted et al. (2022)</td></tr><tr><td>Radius (R<sub>pl</sub>)</td><td>3.85</td><td>R<sub>⊕</sub></td><td>Maxted et al. (2022)</td></tr><tr><td>Orbital period (P<sub>orb</sub>)</td><td>≈2.644</td><td>days</td><td>Lanotte et al. (2014)</td></tr><tr><td>Semi major axis (a<sub>orb</sub>)</td><td>14.56</td><td>R<sub>★</sub></td><td>Maxted et al. (2022)</td></tr><tr><td>...</td><td>≈0.028</td><td>au</td><td>...</td></tr><tr><td>Spin-orbit angle (Ψ)</td><td>103.2</td><td>degrees</td><td>Bourrier et al. (2022)</td></tr></tbody></table>

(Mauas et al. 1997; Peacock et al. 2019), which would imply transition region densities of a few times $10^9$ cm$^{-3}$ at a temperature of $5 \times 10^4$ K. Depending on which chromospheric model we use to guide our assumption of $n_{tr}$, we could have transition region densities that vary from ~$10^9$ to $10^{11}$ cm$^{-3}$. O Fionnagain et al. (2021) explored the effects that changing this parameter had on the wind of a sub-giant star $\lambda$ And and they found that varying $n_{tr}$ by a factor of 100 led to similar wind properties (compare their models A2 and B0, or B1 and D0, from their Table 2). We therefore do not expect that assuming a different value for $n_{tr}$ within the quoted range would change our wind solutions substantially.

As they propagate through magnetic field lines, Alfvén waves transfer their energy and momentum to the wind plasma, thus causing its heating and acceleration. In AWSOM, the wave dissipation is modelled through a turbulent cascade process. Its correlation length is assumed to increase with flux tube radius; thus $L_{\perp} = L_{\perp,0} \sqrt{B_{\star}/B}$. Here, we assume the proportionality factor is the same as that used in the AWSOM solar wind models, namely $L_{\perp,0} \sqrt{B_{\star}} = 1.5 \times 10^5$ m $\sqrt{T} = 1.5 \times 10^4$ km $\sqrt{G}$. Based on the parametric study of O Fionnagain et al. (2021), we also do not expect that changing the value of this base quantity by a factor of ~10 will substantially affect the general properties of our wind solution (see also Shoda et al. 2018, who investigated the effects of the correlation length on wind properties such as its maximum temperature, mass-loss rate, and velocity).

An important quantity of our model is the boundary value of the Alfvén wave flux:

$$S_{A,0} = \rho_{\text{tr}} \langle v_{\perp,0}^2 \rangle v_{A,0}, \qquad (1)$$

where $\langle v_{\perp,0}^2 \rangle^{1/2}$ is the average velocity perturbation of the Alfvén wave, $\rho_{tr} = m_p n_{tr}$ is the mass density at the transition region, and the Alfvén velocity is $v_{A,0} = B_{\star} / \sqrt{4\pi\rho_{tr}}$, where the index '0' indicates value computed at the boundary of the simulation. Here, $B_{\star}$ is the local value (i.e. at each longitude and latitude) of the surface magnetic field as derived in the surface magnetic map (Fig. 1). The previous equation can be rewritten as

$$\frac{\langle v_{\perp,0}^2 \rangle^{1/2}}{65 \text{ km s}^{-1}} = \left( \frac{5 \times 10^{11}}{n_{\text{tr}}[\text{cm}^{-3}]} \right)^{1/4} \left( \frac{S_{A,0}/B_{\star}[\text{erg/(cm}^2 \text{ s G})]}{1.1 \times 10^7} \right)^{1/2} . \quad (2)$$

We note that the ratio $S_{A,0}/B_{\star}$ is assumed to be uniform at the inner boundary, implying that the average velocity amplitude of the Alfvén wave and the energy density of the wave ($\rho_{tr}\langle v_{\perp,0}^2 \rangle$)

are also uniform at the inner boundary. The choice of solar wave flux adopted in previous models of the solar wind with AWSOM (Oran et al. 2013, 2017; van der Holst et al. 2014) was guided by Hinode observations of velocity fluctuations at the chromospheric region (De Pontieu et al. 2007). Indeed, if one assumes that the non-thermal motions of ions at the transition region are associated with Alfvén waves (e.g. Oran et al. 2017), then one way to constrain the Alfvén wave properties in other stars is by measuring the non-thermal broadening of lines formed in the transition region (Boro Saikia et al. 2023). In this case, the non-thermal velocity $\xi$ is related to the average velocity perturbation of the Alfvén wave as $\xi^2 = \frac{1}{2}\langle v_{\perp,0}^2 \rangle$.

The ratio $S_{A,0}/B_{\star}$ is the most important free parameter in our model, and its value affects the wind temperature, its mass-loss rate, and the size of the Alfvén surface, which in turn affects the angular momentum loss of the star (Boro Saikia et al. 2020; O Fionnagain et al. 2021; Kavanagh et al. 2021). Given that, we ran two wind models with two different values of $S_{A,0}/B_{\star}$: $1.1 \times 10^6$ erg (cm$^2$ s G)$^{-1}$ and $1.1 \times 10^7$ erg (cm$^2$ s G)$^{-1}$. The latter is the same value as used in AWSOM solar wind simulations (Sokolov et al. 2013; van der Holst et al. 2014). From now on, we refer to these models as 'I' and 'II'. As we show in Sect. 3, our derived mass-loss rates and computed non-thermal velocities of Model I show better agreement with observations, indicating that the adopted lower Alfvén wave flux is more appropriate to describe the wind of GJ 436.

For GJ 436, which has a surface radial field reaching about 28 G, our choice of parameters translates to correlation lengths of $L_{\perp,0} \approx 2800$ km$^1$ and surface Alfvén fluxes of up to $S_{A,0} \approx [3 \times 10^7, 3 \times 10^8]$ erg cm$^{-2}$s$^{-1}$ for Models I and II, respectively. The surface amplitudes of the wave perturbation are $21$ km s$^{-1}$ (Model I) and $65$ km s$^{-1}$ (Model II), uniform at the lower boundary. The corresponding surface energy densities of the

<footnote n="1">¹ The scaling of the correlation length $L_{\perp,0} \sqrt{B_{\star}}$ has been related to the distance between magnetic flux tubes on the solar surface (Hollweg 1986). Hollweg (1986) empirically estimated a value of $L_{\perp,0} \sqrt{B_{\star}} = 7520$ km $\sqrt{G}$ for the Sun, or $L_{\perp,0} \sim 4300$ km for ~3 G field. We do not know how $L_{\perp,0}$ would change for GJ 436, but because of its smaller radius, we naively expect that the distance between magnetic flux tubes would be smaller than solar for GJ 436, consistent to the value of ~2800 km we use. Based on the parametric study of O Fionnagain et al. (2021), we do not expect a variation in this parameter by a factor of ~10 would substantially affect the general properties of our wind solution.</footnote>
<table><caption>**Table 1.** Parameters adopted in our models for GJ 436 and GJ 436b.</caption><thead><tr><th>Parameter (Star)</th><th>Value</th><th>Units</th><th>Reference</th></tr></thead><tbody><tr><td>Mass (M<sub>★</sub>)</td><td>≈0.45</td><td>M<sub>☉</sub></td><td>Maxted et al. (2022)</td></tr><tr><td>Radius (R<sub>★</sub>)</td><td>≈0.42</td><td>R<sub>☉</sub></td><td>Maxted et al. (2022)</td></tr><tr><td>Rotation period (P<sub>rot</sub>)</td><td>≈44.1</td><td>days</td><td>Bourrier et al. (2018)</td></tr><tr><td>Stellar magnetism</td><td>cf. Fig. 1</td><td></td><td>Bellotti et al. (2023)</td></tr><tr><td>Age</td><td>6<sup>+4</sup><sub>-5</sub></td><td>Gyr</td><td>Torres (2007)</td></tr></tbody></table><table><thead><tr><th>Parameter (Planet)</th><th>Value</th><th>Units</th><th>Reference</th></tr></thead><tbody><tr><td>Mass (M<sub>pl</sub>)</td><td>21.72</td><td>M<sub>⊕</sub></td><td>Maxted et al. (2022)</td></tr><tr><td>Radius (R<sub>pl</sub>)</td><td>3.85</td><td>R<sub>⊕</sub></td><td>Maxted et al. (2022)</td></tr><tr><td>Orbital period (P<sub>orb</sub>)</td><td>≈2.644</td><td>days</td><td>Lanotte et al. (2014)</td></tr><tr><td>Semi major axis (a<sub>orb</sub>)</td><td>14.56</td><td>R<sub>★</sub></td><td>Maxted et al. (2022)</td></tr><tr><td>...</td><td>≈0.028</td><td>au</td><td>...</td></tr><tr><td>Spin-orbit angle (Ψ)</td><td>103.2</td><td>degrees</td><td>Bourrier et al. (2022)</td></tr></tbody></table>

(Mauas et al. 1997; Peacock et al. 2019), which would imply transition region densities of a few times $10^9$ cm$^{-3}$ at a temperature of $5 \times 10^4$ K. Depending on which chromospheric model we use to guide our assumption of $n_{tr}$, we could have transition region densities that vary from ~$10^9$ to $10^{11}$ cm$^{-3}$. O Fionnagain et al. (2021) explored the effects that changing this parameter had on the wind of a sub-giant star $\lambda$ And and they found that varying $n_{tr}$ by a factor of 100 led to similar wind properties (compare their models A2 and B0, or B1 and D0, from their Table 2). We therefore do not expect that assuming a different value for $n_{tr}$ within the quoted range would change our wind solutions substantially.

As they propagate through magnetic field lines, Alfvén waves transfer their energy and momentum to the wind plasma, thus causing its heating and acceleration. In AWSOM, the wave dissipation is modelled through a turbulent cascade process. Its correlation length is assumed to increase with flux tube radius; thus $L_{\perp} = L_{\perp,0} \sqrt{B_{\star}/B}$. Here, we assume the proportionality factor is the same as that used in the AWSOM solar wind models, namely $L_{\perp,0} \sqrt{B_{\star}} = 1.5 \times 10^5$ m $\sqrt{T} = 1.5 \times 10^4$ km $\sqrt{G}$. Based on the parametric study of O Fionnagain et al. (2021), we also do not expect that changing the value of this base quantity by a factor of ~10 will substantially affect the general properties of our wind solution (see also Shoda et al. 2018, who investigated the effects of the correlation length on wind properties such as its maximum temperature, mass-loss rate, and velocity).

An important quantity of our model is the boundary value of the Alfvén wave flux:

$$S_{A,0} = \rho_{\text{tr}} \langle v_{\perp,0}^2 \rangle v_{A,0}, \qquad (1)$$

where $\langle v_{\perp,0}^2 \rangle^{1/2}$ is the average velocity perturbation of the Alfvén wave, $\rho_{tr} = m_p n_{tr}$ is the mass density at the transition region, and the Alfvén velocity is $v_{A,0} = B_{\star} / \sqrt{4\pi\rho_{tr}}$, where the index '0' indicates value computed at the boundary of the simulation. Here, $B_{\star}$ is the local value (i.e. at each longitude and latitude) of the surface magnetic field as derived in the surface magnetic map (Fig. 1). The previous equation can be rewritten as

$$\frac{\langle v_{\perp,0}^2 \rangle^{1/2}}{65 \text{ km s}^{-1}} = \left( \frac{5 \times 10^{11}}{n_{\text{tr}}[\text{cm}^{-3}]} \right)^{1/4} \left( \frac{S_{A,0}/B_{\star}[\text{erg/(cm}^2 \text{ s G})]}{1.1 \times 10^7} \right)^{1/2} . \quad (2)$$

We note that the ratio $S_{A,0}/B_{\star}$ is assumed to be uniform at the inner boundary, implying that the average velocity amplitude of the Alfvén wave and the energy density of the wave ($\rho_{tr}\langle v_{\perp,0}^2 \rangle$)

are also uniform at the inner boundary. The choice of solar wave flux adopted in previous models of the solar wind with AWSOM (Oran et al. 2013, 2017; van der Holst et al. 2014) was guided by Hinode observations of velocity fluctuations at the chromospheric region (De Pontieu et al. 2007). Indeed, if one assumes that the non-thermal motions of ions at the transition region are associated with Alfvén waves (e.g. Oran et al. 2017), then one way to constrain the Alfvén wave properties in other stars is by measuring the non-thermal broadening of lines formed in the transition region (Boro Saikia et al. 2023). In this case, the non-thermal velocity $\xi$ is related to the average velocity perturbation of the Alfvén wave as $\xi^2 = \frac{1}{2}\langle v_{\perp,0}^2 \rangle$.

The ratio $S_{A,0}/B_{\star}$ is the most important free parameter in our model, and its value affects the wind temperature, its mass-loss rate, and the size of the Alfvén surface, which in turn affects the angular momentum loss of the star (Boro Saikia et al. 2020; O Fionnagain et al. 2021; Kavanagh et al. 2021). Given that, we ran two wind models with two different values of $S_{A,0}/B_{\star}$: $1.1 \times 10^6$ erg (cm$^2$ s G)$^{-1}$ and $1.1 \times 10^7$ erg (cm$^2$ s G)$^{-1}$. The latter is the same value as used in AWSOM solar wind simulations (Sokolov et al. 2013; van der Holst et al. 2014). From now on, we refer to these models as 'I' and 'II'. As we show in Sect. 3, our derived mass-loss rates and computed non-thermal velocities of Model I show better agreement with observations, indicating that the adopted lower Alfvén wave flux is more appropriate to describe the wind of GJ 436.

For GJ 436, which has a surface radial field reaching about 28 G, our choice of parameters translates to correlation lengths of $L_{\perp,0} \approx 2800$ km$^1$ and surface Alfvén fluxes of up to $S_{A,0} \approx [3 \times 10^7, 3 \times 10^8]$ erg cm$^{-2}$s$^{-1}$ for Models I and II, respectively. The surface amplitudes of the wave perturbation are $21$ km s$^{-1}$ (Model I) and $65$ km s$^{-1}$ (Model II), uniform at the lower boundary. The corresponding surface energy densities of the

<footnote n="1">¹ The scaling of the correlation length $L_{\perp,0} \sqrt{B_{\star}}$ has been related to the distance between magnetic flux tubes on the solar surface (Hollweg 1986). Hollweg (1986) empirically estimated a value of $L_{\perp,0} \sqrt{B_{\star}} = 7520$ km $\sqrt{G}$ for the Sun, or $L_{\perp,0} \sim 4300$ km for ~3 G field. We do not know how $L_{\perp,0}$ would change for GJ 436, but because of its smaller radius, we naively expect that the distance between magnetic flux tubes would be smaller than solar for GJ 436, consistent to the value of ~2800 km we use. Based on the parametric study of O Fionnagain et al. (2021), we do not expect a variation in this parameter by a factor of ~10 would substantially affect the general properties of our wind solution.</footnote>
<figure><img src="image_2.png" /><figcaption>**Fig. 2.** 3D view of the magnetised stellar wind for Models I and II. The surface magnetic field from March-June 2016 reconstructed by Bellotti et al. (2023) is used at the boundary of the models. The grey streamlines are the magnetic field lines embedded in the wind. The rotation axis of the star is along the $z_\star$-axis. The 2D cut is placed at the $x_\star z_\star$ plane, and for reference the orbital distance of the planet is shown with the translucent sphere at $14.56 R_\star \approx 0.028$ au. The contour shows the wind velocity, and the Alfvén surface is shown by the translucent outer surface. Models I and II predict stellar wind mass-loss rates of $1.1 \times 10^{-15} M_\odot \text{ yr}^{-1}$ and $2.5 \times 10^{-14} M_\odot \text{ yr}^{-1}$.</figcaption></figure>

wave are 11 and 35 erg cm$^{-3}$. With these parameters, the predicted non-thermal velocities at the transition region are $\xi \approx 14$ and $45 \text{ km s}^{-1}$ for Models I and II, respectively. The other remaining parameters of the model, namely the stochastic heating parameter, the collisionless heat conduction parameter and the collisionless heat flux parameter have the same values as in Table 1 of van der Holst et al. (2014). The energy equations are solved in their non-conservative form (i.e. solving for the pressure instead of total energy).

Figure 2 shows the 3D outputs of our simulations, after they have reached steady state. The grey streamlines represent the stellar wind magnetic field lines; this geometry results from the self-consistent interaction between field lines and the stellar wind flow. The rotation axis lies along positive $z_\star$ – we note that the magnetic field geometry is approximately that of an aligned dipole close to the stellar surface, with elongated closed field lines at large distances and open fields around the stellar rotation poles. The stellar wind mass-loss rate we obtained are $1.1 \times 10^{-15} M_\odot \text{ yr}^{-1}$ and $2.5 \times 10^{-14} M_\odot \text{ yr}^{-1}$ for Models I and II, respectively. For comparison, the solar wind mass-loss rate is $2 \times 10^{-14} M_\odot \text{ yr}^{-1}$ – given that all of our input values for the Alfvén wave and transition region are solar for Model II, and the surface magnetic flux of GJ 436 is similar to solar<fnref n="2" />, it is not surprising that the mass-loss rate derived in Model II is similar to that of the solar wind. In Sect. 3 we perform a more detailed comparison of our derived mass-loss rates with those derived in other studies of GJ 436. Figure 2 also shows a cut in

the $y_\star z_\star$ plane, where colours indicate the stellar wind speed, where we see that along open field lines the wind is accelerated more quickly, with low-speed streams around closed field line regions. As expected, higher velocities are achieved in Model II, with higher Alfvén wave fluxes. The Alfvén surface is shown by the translucent surface. For context, we plot a white circumference at the orbital distance of GJ 436b (for simplicity, we assume a circular orbit). We note most of the time the planet would be in the region of open magnetic field lines and that the orbital distance is nearly always within the Alfvén surface for both models, which can have implications for star-planet interactions (Folsom et al. 2020) and induced radio emission on the star (Kavanagh et al. 2021, 2022), for example. We come back to this in Sect. 4.

### 3. Space weather conditions of GJ 436b

The space weather of GJ 436b is formed by the particles and ambient magnetic field of the stellar outflow. Here, we compute the properties of stellar wind at the position of GJ 436b. Because GJ 436b has a very misaligned orbit, we first compute the position of the planet as a function of time, following the coordinate transformations presented in Appendix A. To calculate the phase of the orbit with respect to the rotational phase of the star, we consider the transit midpoint at $HJD_t = 2458947.26212(12)$ (Maxted et al. 2022), with an orbital period of 2.64389803(27) d (Lanotte et al. 2014). At stellar rotation phase zero ($HJD_s = 2457464.4967$, corresponding to the first spectropolarimetric observation (Bellotti et al. 2023), and $\{x_\star, y_\star\} = \{R_\star, 0\}$ in our simulation coordinates), the planet was at an orbital phase of $\varphi_0 = 0.1745$ after the preceding transit midpoint. GJ 436b has a nearly polar orbit with its orbital axis and stellar spin misaligned by $\Psi = 103.2^\circ$ (Bourrier et al. 2022). Figure 3 shows the trajectories of the planet as seen in the reference frame of the star during one rotational period of the star (about 44 days), during which the planet completes about 16 orbital revolutions.

<footnote n="2">² From the magnetic map of GJ 436, we find a surface magnetic flux of $\Phi_0 = 1.4 \times 10^{23}$ Mx, which is very similar to solar values (for an average large-scale field strength of ~3 G, the solar surface flux is $\Phi_0 = 3 \times 4\pi R_\star^2 \approx 1.8 \times 10^{23}$ Mx). In stellar wind theory, the mass-loss rate is linearly proportional to the open magnetic flux ($\Phi_{\text{open}}$), which is also proportional to $\Phi_0$ (e.g. Vidotto et al. 2014). This is also seen in AWSOM models of winds of Sun-like stars (Evensberget et al. 2023). Because Model II uses the same solar wind parameters and GJ 436 has a similar surface flux to the Sun, this model produces a mass-loss rate that is thus similar to that of the Sun.</footnote>
<figure><img src="image_3.png" /><figcaption>**Fig. 2.** 3D view of the magnetised stellar wind for Models I and II. The surface magnetic field from March-June 2016 reconstructed by Bellotti et al. (2023) is used at the boundary of the models. The grey streamlines are the magnetic field lines embedded in the wind. The rotation axis of the star is along the $z_{\star}$-axis. The 2D cut is placed at the $x_{\star}z_{\star}$ plane, and for reference the orbital distance of the planet is shown with the translucent sphere at $14.56 R_{\star} \approx 0.028$ au. The contour shows the wind velocity, and the Alfvén surface is shown by the translucent outer surface. Models I and II predict stellar wind mass-loss rates of $1.1 \times 10^{-15} M_{\odot} \text{ yr}^{-1}$ and $2.5 \times 10^{-14} M_{\odot} \text{ yr}^{-1}$.</figcaption></figure>

wave are 11 and 35 erg cm$^{-3}$. With these parameters, the predicted non-thermal velocities at the transition region are $\xi \approx 14$ and $45 \text{ km s}^{-1}$ for Models I and II, respectively. The other remaining parameters of the model, namely the stochastic heating parameter, the collisionless heat conduction parameter and the collisionless heat flux parameter have the same values as in Table 1 of van der Holst et al. (2014). The energy equations are solved in their non-conservative form (i.e. solving for the pressure instead of total energy).

Figure 2 shows the 3D outputs of our simulations, after they have reached steady state. The grey streamlines represent the stellar wind magnetic field lines; this geometry results from the self-consistent interaction between field lines and the stellar wind flow. The rotation axis lies along positive $z_{\star}$ – we note that the magnetic field geometry is approximately that of an aligned dipole close to the stellar surface, with elongated closed field lines at large distances and open fields around the stellar rotation poles. The stellar wind mass-loss rate we obtained are $1.1 \times 10^{-15} M_{\odot} \text{ yr}^{-1}$ and $2.5 \times 10^{-14} M_{\odot} \text{ yr}^{-1}$ for Models I and II, respectively. For comparison, the solar wind mass-loss rate is $2 \times 10^{-14} M_{\odot} \text{ yr}^{-1}$ – given that all of our input values for the Alfvén wave and transition region are solar for Model II, and the surface magnetic flux of GJ 436 is similar to solar<fnref n="2" />, it is not surprising that the mass-loss rate derived in Model II is similar to that of the solar wind. In Sect. 3 we perform a more detailed comparison of our derived mass-loss rates with those derived in other studies of GJ 436. Figure 2 also shows a cut in

the $y_{\star}z_{\star}$ plane, where colours indicate the stellar wind speed, where we see that along open field lines the wind is accelerated more quickly, with low-speed streams around closed field line regions. As expected, higher velocities are achieved in Model II, with higher Alfvén wave fluxes. The Alfvén surface is shown by the translucent surface. For context, we plot a white circumference at the orbital distance of GJ 436b (for simplicity, we assume a circular orbit). We note most of the time the planet would be in the region of open magnetic field lines and that the orbital distance is nearly always within the Alfvén surface for both models, which can have implications for star-planet interactions (Folsom et al. 2020) and induced radio emission on the star (Kavanagh et al. 2021, 2022), for example. We come back to this in Sect. 4.

### 3. Space weather conditions of GJ 436b

The space weather of GJ 436b is formed by the particles and ambient magnetic field of the stellar outflow. Here, we compute the properties of stellar wind at the position of GJ 436b. Because GJ 436b has a very misaligned orbit, we first compute the position of the planet as a function of time, following the coordinate transformations presented in Appendix A. To calculate the phase of the orbit with respect to the rotational phase of the star, we consider the transit midpoint at $HJD_t = 2458947.26212(12)$ (Maxted et al. 2022), with an orbital period of 2.64389803(27) d (Lanotte et al. 2014). At stellar rotation phase zero ($HJD_s = 2457464.4967$, corresponding to the first spectropolarimetric observation (Bellotti et al. 2023), and $\{x_{\star}, y_{\star}\} = \{R_{\star}, 0\}$ in our simulation coordinates), the planet was at an orbital phase of $\varphi_0 = 0.1745$ after the preceding transit midpoint. GJ 436b has a nearly polar orbit with its orbital axis and stellar spin misaligned by $\Psi = 103.2^\circ$ (Bourrier et al. 2022). Figure 3 shows the trajectories of the planet as seen in the reference frame of the star during one rotational period of the star (about 44 days), during which the planet completes about 16 orbital revolutions.

<footnote n="2">² From the magnetic map of GJ 436, we find a surface magnetic flux of $\Phi_0 = 1.4 \times 10^{23}$ Mx, which is very similar to solar values (for an average large-scale field strength of ~3 G, the solar surface flux is $\Phi_0 = 3 \times 4\pi R_{\star}^2 \approx 1.8 \times 10^{23}$ Mx). In stellar wind theory, the mass-loss rate is linearly proportional to the open magnetic flux ($\Phi_{\text{open}}$), which is also proportional to $\Phi_0$ (e.g. Vidotto et al. 2014). This is also seen in AWSOM models of winds of Sun-like stars (Evensberget et al. 2023). Because Model II uses the same solar wind parameters and GJ436 has a similar surface flux to the Sun, this model produces a mass-loss rate that is thus similar to that of the Sun.</footnote>
<figure><img src="image_3.png" /><figcaption>**Fig. 3.** Trajectories of the planet as seen in the reference frame of the star during one rotational period of the star (about 44 days; see the colour bar), during which the planet completes about 16 orbital revolutions. The grey circles indicate the approximate times of mid-transit.</figcaption></figure>

As presented in Appendix A, in the case of GJ 436b, which has a highly misaligned orbit, the time it takes for the planet to be at the same sub-stellar point is a complex combination of orbital and rotation periods, as well as the obliquity of the orbit. Because of this, the local condition of the stellar wind is constantly changing. Figure 4 shows the stellar wind proton density, velocity, total pressure (ram + thermal + magnetic) and magnetic field strength, all computed at the orbit of GJ 436b. The solid lines are for Model I, and the dashed lines for Model II. Colour represents the time evolution following the same colour scheme as shown in Fig. 3. We note that to compute the ram pressure, we use the relative velocity between the stellar wind velocity **u** and the Keplerian motion of the planet **v**<sub>K</sub>: Δ**u** = **u** − **v**<sub>K</sub>. For GJ 436b, the Keplerian velocity is **v**<sub>K</sub> = 118(sin Ψ̂θ + cos Ψ̂φ) km s<sup>−1</sup>, where we adopted a stellar mass of ≈0.45 M<sub>☉</sub> and an orbital distance of 14.56 R<sub>⋆</sub> ≈ 0.028 au (Maxted et al. 2022).

The range of wind densities at the orbit of the planet is
[0.12, 5.0] × 10³ cm⁻³ for Model I and [2.8, 81] × 10³ cm⁻³ for
Model II. Their time averages over one stellar rotation period are
510 and 6200 cm⁻³ for Models I and II, respectively (horizon-
tal grey lines shown in Fig. 4). Similarly, for wind velocities,
we found that the range of values that the planet experience
are [30, 770] km s⁻¹ and [220, 830] km s⁻¹ for Models I and
II, respectively, with averages of 510 and 640 km s⁻¹. When
modelling the transmission spectroscopic transit of GJ 436b,
Bourrier et al. (2016) derived densities of 2<sup>+2.2</sup><sub>-1.2</sub> × 10³ cm⁻³ and
wind speeds of 85<sup>+6</sup><sub>-16</sub> km s⁻¹. These observations are taken dur-
ing transit, so the values they derived should be compared to our
model values around mid-transit times. For this comparison, we
added to the top panels of Fig. 4 crosses at the approximate mid-
transit times at the same level of ion density (2 × 10³ cm⁻³) and
speed (85 km s⁻¹) derived in the work of Bourrier et al. (2016).
We see that these crosses overlap really well with values of
Model I (solid lines), but not with Model II (dashed lines), imply-
ing that Model I, with the lowest surface flux of Alfvén waves,
better matches the spectroscopic transit values from Bourrier
et al. (2016). Indeed, from these observations, the estimated

mass-loss rate of 1.2<sup>+1.3</sup><sub>-0.75</sub> × 10<sup>-15</sup> M<sub>☉</sub> yr<sup>-1</sup> (Vidotto & Bourrier 2017) better agree with the values derived in Model I (1.1 × 10<sup>-15</sup> M<sub>☉</sub> yr<sup>-1</sup>), but not with Model II (2.5 × 10<sup>-14</sup> M<sub>☉</sub> yr<sup>-1</sup>). Our chosen values of S<sub>A,0</sub>/B<sub>⋆</sub> correspond to non-thermal velocities ξ ≈ 14 and 45 km s<sup>-1</sup> at the transition region (see Sect. 2). Recently, Boro Saikia et al. (2023) reported ξ ranging from ~5 to ~30 km s<sup>-1</sup> for stars with similar rotation periods and effective temperatures as GJ 436. In this case too, the non-thermal velocities from Model I have a better agreement with the range reported in Boro Saikia et al. (2023). Altogether, our derived mass-loss rates and non-thermal velocities suggest that the ‘canonical’ value of S<sub>A,0</sub>/B<sub>⋆</sub> = 1.1 × 10<sup>7</sup> erg (cm<sup>2</sup> s G)<sup>-1</sup>, derived from solar wind models, might not be applicable to the winds of (old) M dwarf stars.

The bottom-right panel of Fig. 4 shows the stellar wind magnetic field strength at the orbit of GJ 436b. We find an average stellar magnetic field strength of 0.013 G and 0.026 G for Models I and II, respectively. We use these values to then compute the magnetic pressure (B²/8π). We also compute the ram pressure of the stellar wind (ρΔu²). Adding these two pressures to the thermal pressure, we obtain the total pressure of the stellar wind at the position of the planet, which is shown at the bottom-left panel of Fig. 4. We find average total pressures of 8.2 × 10⁻⁶ dyn/cm⁻² and 59 × 10⁻⁶ dyn/cm⁻² for Models I and II, respectively. These values will be used in Sect. 4, when we compute the strength of star-planet interactions and the size of the magnetosphere of the planet.

**4. Stellar wind–planet interactions**

Because of its close proximity to the star, and its sub-Alfvénic
motion, GJ 436b can magnetically interact with its host star.
Recently, Loyd et al. (2023) studied the presence of star–planet
interaction signatures in GJ 436 in ultraviolet spectroscopic
observations. While they did not find line fluxes modulated with
the orbital period of GJ 436b (traditionally used as evidence
of star–planet interactions; e.g. Pineda & Villadsen 2023), they
found an enhancement of low energy flares and a lack of more
energetic flares in GJ 436, compared to other M dwarfs. They
suggested that this could be caused by the presence of GJ 436b,
whereby the reduction of the number of high energy flares in
GJ 436 would be caused by the planet triggering early flaring
events (i.e. before enough energy would have been built up in
a larger stellar flare). As a consequence, high-energy flares in
GJ 436 would not have time to store energy and would be pre-
maturely released, explaining the excess of low energy flares and
lack of high energy flares in the system. In their analysis, they
also derived a maximum power for the star–planet interaction,
detected in the far ultraviolet lines, of ~3 × 10²⁴ erg s⁻¹. Convert-
ing this number to a bolometric power, they found an upper limit
of ≤10²⁶ erg s⁻¹ for the maximum power released in star–planet
interactions in the GJ 436 system.

To compare with the derived power from Loyd et al. (2023),
in this section we compute the strength of star–planet interac-
tions. Here, we use two different scenarios for this calculation.
One scenario is based on the idea that magnetic reconnection
occurring between stellar and planetary magnetic field lines
releases magnetic energy (Jardine & Cameron 2008; Lanza
2009; Vidotto et al. 2010). This energy can accelerate elec-
trons, which travel towards the star. This idea has been used
to explain anomalous hot spots at stellar chromospheres (Lanza
2009), as reported in observational works (e.g. Shkolnik et al.
2005; Pagano et al. 2009; Cauley et al. 2019). The other
<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Trajectories of the planet as seen in the reference frame of the star during one rotational period of the star (about 44 days; see the colour bar), during which the planet completes about 16 orbital revolutions. The grey circles indicate the approximate times of mid-transit.</figcaption></figure>

As presented in Appendix A, in the case of GJ 436b, which has a highly misaligned orbit, the time it takes for the planet to be at the same sub-stellar point is a complex combination of orbital and rotation periods, as well as the obliquity of the orbit. Because of this, the local condition of the stellar wind is constantly changing. Figure 4 shows the stellar wind proton density, velocity, total pressure (ram + thermal + magnetic) and magnetic field strength, all computed at the orbit of GJ 436b. The solid lines are for Model I, and the dashed lines for Model II. Colour represents the time evolution following the same colour scheme as shown in Fig. 3. We note that to compute the ram pressure, we use the relative velocity between the stellar wind velocity **u** and the Keplerian motion of the planet **v**<sub>K</sub>: $\Delta\mathbf{u} = \mathbf{u} - \mathbf{v}_K$. For GJ 436b, the Keplerian velocity is $\mathbf{v}_K = 118(\sin\Psi\hat{\theta} + \cos\Psi\hat{\phi})$ km s<sup>-1</sup>, where we adopted a stellar mass of ≈0.45 $M_\odot$ and an orbital distance of 14.56 $R_\star$ ≈ 0.028 au (Maxted et al. 2022).

The range of wind densities at the orbit of the planet is
[0.12, 5.0] × 10³ cm⁻³ for Model I and [2.8, 81] × 10³ cm⁻³ for
Model II. Their time averages over one stellar rotation period are
510 and 6200 cm⁻³ for Models I and II, respectively (horizon-
tal grey lines shown in Fig. 4). Similarly, for wind velocities,
we found that the range of values that the planet experience
are [30, 770] km s⁻¹ and [220, 830] km s⁻¹ for Models I and
II, respectively, with averages of 510 and 640 km s⁻¹. When
modelling the transmission spectroscopic transit of GJ 436b,
Bourrier et al. (2016) derived densities of 2<sup>+2.2</sup>/<sub>-1.2</sub> × 10³ cm⁻³ and
wind speeds of 85<sup>+6</sup>/<sub>-16</sub> km s⁻¹. These observations are taken dur-
ing transit, so the values they derived should be compared to our
model values around mid-transit times. For this comparison, we
added to the top panels of Fig. 4 crosses at the approximate mid-
transit times at the same level of ion density (2 × 10³ cm⁻³) and
speed (85 km s⁻¹) derived in the work of Bourrier et al. (2016).
We see that these crosses overlap really well with values of
Model I (solid lines), but not with Model II (dashed lines), imply-
ing that Model I, with the lowest surface flux of Alfvén waves,
better matches the spectroscopic transit values from Bourrier
et al. (2016). Indeed, from these observations, the estimated

mass-loss rate of 1.2<sup>+1.3</sup>/<sub>-0.75</sub> × 10<sup>-15</sup> M<sub>☉</sub> yr<sup>-1</sup> (Vidotto & Bourrier 2017) better agree with the values derived in Model I (1.1 × 10<sup>-15</sup> M<sub>☉</sub> yr<sup>-1</sup>), but not with Model II (2.5 × 10<sup>-14</sup> M<sub>☉</sub> yr<sup>-1</sup>). Our chosen values of S<sub>A,0</sub>/B<sub>⋆</sub> correspond to non-thermal velocities ξ ≈ 14 and 45 km s<sup>-1</sup> at the transition region (see Sect. 2). Recently, Boro Saikia et al. (2023) reported ξ ranging from ~5 to ~30 km s<sup>-1</sup> for stars with similar rotation periods and effective temperatures as GJ 436. In this case too, the non-thermal velocities from Model I have a better agreement with the range reported in Boro Saikia et al. (2023). Altogether, our derived mass-loss rates and non-thermal velocities suggest that the ‘canonical’ value of S<sub>A,0</sub>/B<sub>⋆</sub> = 1.1 × 10<sup>7</sup> erg (cm<sup>2</sup> s G)<sup>-1</sup>, derived from solar wind models, might not be applicable to the winds of (old) M dwarf stars.

The bottom-right panel of Fig. 4 shows the stellar wind magnetic field strength at the orbit of GJ 436b. We find an average stellar magnetic field strength of 0.013 G and 0.026 G for Models I and II, respectively. We use these values to then compute the magnetic pressure ($B^2/8\pi$). We also compute the ram pressure of the stellar wind ($\rho\Delta u^2$). Adding these two pressures to the thermal pressure, we obtain the total pressure of the stellar wind at the position of the planet, which is shown at the bottom-left panel of Fig. 4. We find average total pressures of $8.2 \times 10^{-6}$ dyn/cm$^{-2}$ and $59 \times 10^{-6}$ dyn/cm$^{-2}$ for Models I and II, respectively. These values will be used in Sect. 4, when we compute the strength of star-planet interactions and the size of the magnetosphere of the planet.

**4. Stellar wind–planet interactions**

Because of its close proximity to the star, and its sub-Alfvénic
motion, GJ 436b can magnetically interact with its host star.
Recently, Loyd et al. (2023) studied the presence of star–planet
interaction signatures in GJ 436 in ultraviolet spectroscopic
observations. While they did not find line fluxes modulated with
the orbital period of GJ 436b (traditionally used as evidence
of star–planet interactions; e.g. Pineda & Villadsen 2023), they
found an enhancement of low energy flares and a lack of more
energetic flares in GJ 436, compared to other M dwarfs. They
suggested that this could be caused by the presence of GJ 436b,
whereby the reduction of the number of high energy flares in
GJ 436 would be caused by the planet triggering early flaring
events (i.e. before enough energy would have been built up in
a larger stellar flare). As a consequence, high-energy flares in
GJ 436 would not have time to store energy and would be pre-
maturely released, explaining the excess of low energy flares and
lack of high energy flares in the system. In their analysis, they
also derived a maximum power for the star–planet interaction,
detected in the far ultraviolet lines, of ~3 × 10<sup>24</sup> erg s<sup>-1</sup>. Convert-
ing this number to a bolometric power, they found an upper limit
of ≤10<sup>26</sup> erg s<sup>-1</sup> for the maximum power released in star–planet
interactions in the GJ 436 system.

To compare with the derived power from Loyd et al. (2023),
in this section we compute the strength of star–planet interac-
tions. Here, we use two different scenarios for this calculation.
One scenario is based on the idea that magnetic reconnection
occurring between stellar and planetary magnetic field lines
releases magnetic energy (Jardine & Cameron 2008; Lanza
2009; Vidotto et al. 2010). This energy can accelerate elec-
trons, which travel towards the star. This idea has been used
to explain anomalous hot spots at stellar chromospheres (Lanza
2009), as reported in observational works (e.g. Shkolnik et al.
2005; Pagano et al. 2009; Cauley et al. 2019). The other
<figure><img src="image_4.png" /><figcaption>Fig. 4. Stellar wind proton density, velocity, total pressure (sum of ram, including the planet's own orbital motion and the thermal and magnetic pressures), and total magnetic field strength at the orbit of GJ 436b. The solid lines are for Model I and the dashed lines for Model II. Colour represents the time evolution and follows the colour scheme shown in Fig. 3. The horizontal grey lines are the time average of the stellar wind properties over one stellar rotation period (about 16 orbital periods of the planet). The grey arrows in the upper horizontal axes indicate approximate times of mid-transits. In the upper panels, the crosses indicate the wind values derived in the work of Bourrier et al. (2016), and they are plotted at the mid-transit times. Note how they overlap really well with the values derived in our Model I, indicating that Model I provides a better description of the wind of GJ 436.</figcaption></figure>

scenario we also explore is based on the idea that a planet moving at sub-Alfvénic speed through the magnetised wind of its host star can trigger MHD Alfvén waves, which can travel towards the star along magnetic field lines (Preusse et al. 2006; Saur et al. 2013; Strugarek 2016). As they travel, these waves dissipate their energies, which are radiated away.

We next estimated the power released in the star–planet interaction using the two scenarios mentioned above. To compute the strength of star–planet interactions in both scenarios, we started by calculating the (unsigned) Poynting flux of the stellar wind that impacts on the planet:

$$P = \frac{|(\Delta\mathbf{u} \times \mathbf{B}) \times \mathbf{B}|}{4\pi} = \frac{B^2 \Delta\mathbf{u} \sin\theta}{4\pi} \quad (3)$$

In this expression, $B \sin\theta$ is the component of the stellar wind magnetic field that is perpendicular to $\Delta\mathbf{u}$, $\theta$ is the angle between $\Delta\mathbf{u}$ and $\mathbf{B}$ – all these quantities are computed at the position of the planet, which lies at an orbital distance $a_{\text{orb}} = 14.56 R_{\star}$. We find that the average Poynting fluxes computed at the orbit of GJ 436b are $P \approx 160$ and 650 erg cm$^{-2}$ s$^{-1}$ (or 0.1 to 0.6 W m$^{-2}$)

for Models I and II, respectively. In the context of planets moving at super-Alfvénic speeds, the Poynting flux can be used to estimate the amount of power dissipated in the magnetosphere of a planet, which can then power electron-cyclotron radio emission (e.g. Vidotto et al. 2019). In the case of GJ 436b, though, our wind models indicate that the motion is sub-Alfvénic (cf. Sect. 2) and we use the Poynting flux to compute the energy released in the Alfvén wing scenario and the reconnection scenario. We note that the Poynting flux only depends on the characteristics of the wind at a certain orbital distance and does not depend on the properties of the planet, such as its size or magnetic field strength, contrary to the power released in star–planet interactions, as we see below.

## 4.1. Scenario 1: Alfvén wings

A planet moving at sub-Alfvénic speed through the magnetised wind of its host star can trigger MHD Alfvén waves, which can travel towards the star along magnetic field lines. As they travel, these waves can then dissipate their energies and the maximum
<figure><img src="image_5.png" /><figcaption>Fig. 4. Stellar wind proton density, velocity, total pressure (sum of ram, including the planet's own orbital motion and the thermal and magnetic pressures), and total magnetic field strength at the orbit of GJ 436b. The solid lines are for Model I and the dashed lines for Model II. Colour represents the time evolution and follows the colour scheme shown in Fig. 3. The horizontal grey lines are the time average of the stellar wind properties over one stellar rotation period (about 16 orbital periods of the planet). The grey arrows in the upper horizontal axes indicate approximate times of mid-transits. In the upper panels, the crosses indicate the wind values derived in the work of Bourrier et al. (2016), and they are plotted at the mid-transit times. Note how they overlap really well with the values derived in our Model I, indicating that Model I provides a better description of the wind of GJ 436.</figcaption></figure>

scenario we also explore is based on the idea that a planet moving at sub-Alfvénic speed through the magnetised wind of its host star can trigger MHD Alfvén waves, which can travel towards the star along magnetic field lines (Preusse et al. 2006; Saur et al. 2013; Strugarek 2016). As they travel, these waves dissipate their energies, which are radiated away.

We next estimated the power released in the star–planet interaction using the two scenarios mentioned above. To compute the strength of star–planet interactions in both scenarios, we started by calculating the (unsigned) Poynting flux of the stellar wind that impacts on the planet:

$$P = \frac{|(\Delta\mathbf{u} \times \mathbf{B}) \times \mathbf{B}|}{4\pi} = \frac{B^2 \Delta\mathbf{u} \sin\theta}{4\pi} \quad (3)$$

In this expression, $B \sin\theta$ is the component of the stellar wind magnetic field that is perpendicular to $\Delta\mathbf{u}$, $\theta$ is the angle between $\Delta\mathbf{u}$ and $\mathbf{B}$ – all these quantities are computed at the position of the planet, which lies at an orbital distance $a_{\text{orb}} = 14.56 R_{\star}$. We find that the average Poynting fluxes computed at the orbit of GJ 436b are $P \approx 160$ and 650 erg cm$^{-2}$ s$^{-1}$ (or 0.1 to 0.6 W m$^{-2}$)

for Models I and II, respectively. In the context of planets moving at super-Alfvénic speeds, the Poynting flux can be used to estimate the amount of power dissipated in the magnetosphere of a planet, which can then power electron-cyclotron radio emission (e.g. Vidotto et al. 2019). In the case of GJ 436b, though, our wind models indicate that the motion is sub-Alfvénic (cf. Sect. 2) and we use the Poynting flux to compute the energy released in the Alfvén wing scenario and the reconnection scenario. We note that the Poynting flux only depends on the characteristics of the wind at a certain orbital distance and does not depend on the properties of the planet, such as its size or magnetic field strength, contrary to the power released in star–planet interactions, as we see below.

## 4.1. Scenario 1: Alfvén wings

A planet moving at sub-Alfvénic speed through the magnetised wind of its host star can trigger MHD Alfvén waves, which can travel towards the star along magnetic field lines. As they travel, these waves can then dissipate their energies and the maximum
power that is radiated away can be estimated by integrating the Poynting flux of the stellar wind over the cross-section of two Alfvén 'wings'. In the limit of sub-Alfvénic speeds, with Alfvén-Mach number $M_A = \Delta u/v_A < 1$, the dissipated power is (Saur et al. 2013)

$$ \mathcal{P}_{\text{wg}} \approx P(2\pi r_M^2)\alpha^2 M_A \sin\theta = \frac{\Delta u^2 B \sin^2\theta}{2} r_M^2 \alpha^2 \sqrt{4\pi\rho}. \quad (4) $$

Here, $\alpha$ is a dimensionless parameter related to the conductive properties of the planet and it can be understood as the strength of the interaction; here, we assume $\alpha = 1$, thus deriving the maximum power $\mathcal{P}_{\text{wg}}$. In the expression above, the planet is assumed to be magnetised, with a magnetospheric size $r_M$. Through pressure balance between the magnetic pressure of the planet and the total pressure of the stellar wind (i.e. the sum of the magnetic, thermal, and ram pressures; see the bottom-left panel of Fig. 4), we have (e.g. Vidotto et al. 2014)

$$ \frac{r_M}{R_p} \approx f \left[ \frac{(B_p/2)^2/(8\pi)}{p_{\text{tot}}(a_{\text{orb}})} \right]^{1/6}, \quad (5) $$

where $B_p$ is the polar field strength of the planet's dipolar magnetic field (i.e. twice the value of the equatorial field strength) and $f \approx 2^{2/6}$ is a correction factor used to account for the effects of currents (e.g. Cravens 2004). The value of $B_p$ for GJ 436b is currently poorly constrained: Bourrier et al. (2016) derived an upper limit for the equatorial field strength of ~3 G, consistent with the upper limit estimates of ~10 G from Loyd et al. (2023). Assuming $B_p = 2$ G leads to an average magnetosphere size of 5.2 $R_p$ and 3.7 $R_p$ for stellar wind Models I and II, respectively, as seen on the top panel of Fig. 5. The maximum power that is radiated by Alfvén waves travelling towards the star is shown on the middle panel of this figure, where we see on average that the radiated power is $\mathcal{P}_{\text{wg}} \approx 1.2 \times 10^{22}$ erg s$^{-1}$ and $4.8 \times 10^{22}$ erg s$^{-1}$ for Models I and II, respectively. It is worth noting that the radiated power depends on the assumed magnetic field strength of the planet indirectly through $r_M$. As $r_M \propto B_p^{1/3}$, we have that $\mathcal{P}_{\text{wg}} \propto B_p^{2/3}$. Thus, for a 20 G polar field strength (i.e. 10 G at the equator), the radiated power increases by $10^{2/3}$, which is nearly a factor of 5, resulting in $\mathcal{P}_{\text{wg}} \approx 5.6 \times 10^{22}$ erg s$^{-1}$ and $22 \times 10^{22}$ erg s$^{-1}$ for Models I and II, respectively.

## 4.2. Scenario 2: Star-planet magnetic reconnection

Magnetic reconnection occurs when magnetic field lines of different polarities interact with each other. Magnetic reconnection has been considered in the case of star-planet interactions by several authors (e.g. Jardine & Cameron 2008; Lanza 2009; Vidotto et al. 2010). While Jardine & Cameron (2008); Vidotto et al. (2010) considered reconnection events to power radio emission from exoplanets, Lanza (2009) studied how reconnection events could trigger chromospheric hot spots. We can estimate the power released in the reconnection by integrating the Poynting flux, $P$ (Eq. (3)) over the cross-section of the planet's magnetosphere (Lanza 2009):

$$ \mathcal{P}_{\text{rec}} = \gamma_{\text{rec}} P(\pi R_M^2) = \frac{1}{4} \gamma_{\text{rec}} r_M^2 B^2 \Delta u \sin\theta, \quad (6) $$

where $\gamma_{\text{rec}}$ is an efficiency of the reconnection related to the angle between the stellar and planetary field lines. Even though in the GJ 436 system the planet will move through different orientations of the stellar wind magnetic field (e.g. going through

<figure><img src="image_5.png" /><figcaption>Fig. 5. Strength of the stellar wind-planet interactions computed in this work. Top: magnetospheric size of the planet, assuming a planetary magnetic field strength of 2 G at the pole (i.e., 1 G at the equator). Middle: maximum dissipated power of Alfvén waves travelling towards the star (Saur et al. 2013). Bottom: maximum dissipated power from the star-planet magnetic reconnection scenario (Lanza 2009). The solid lines are for Model I, and the dashed lines for Model II. The colour scheme follows the one adopted in Fig. 3. The grey arrows in the upper horizontal axes indicate approximate times of mid-transits.</figcaption></figure>

closed and open field lines), and thus $\gamma_{\text{rec}}$ would differ at different points in the orbit, for simplicity here we take its maximum value of $\gamma_{\text{rec}} = 1$, so that our power estimates $\mathcal{P}_{\text{rec}}$ are upper limits. The results of this scenario are shown in the bottom panel of Fig. 5, where we see that the average powers released in the reconnection scenario are $\mathcal{P}_{\text{rec}} \approx 8 \times 10^{22}$ erg s$^{-1}$ and $17 \times 10^{22}$ erg s$^{-1}$ for Models I and II, respectively, where we
power that is radiated away can be estimated by integrating the Poynting flux of the stellar wind over the cross-section of two Alfvén 'wings'. In the limit of sub-Alfvénic speeds, with Alfvén-Mach number $M_A = \Delta u/v_A < 1$, the dissipated power is (Saur et al. 2013)

$$ \mathcal{P}_{\text{wg}} \approx P(2\pi r_M^2)\alpha^2 M_A \sin\theta = \frac{\Delta u^2 B \sin^2\theta}{2} r_M^2 \alpha^2 \sqrt{4\pi\rho}. \quad (4) $$

Here, $\alpha$ is a dimensionless parameter related to the conductive properties of the planet and it can be understood as the strength of the interaction; here, we assume $\alpha = 1$, thus deriving the maximum power $\mathcal{P}_{\text{wg}}$. In the expression above, the planet is assumed to be magnetised, with a magnetospheric size $r_M$. Through pressure balance between the magnetic pressure of the planet and the total pressure of the stellar wind (i.e. the sum of the magnetic, thermal, and ram pressures; see the bottom-left panel of Fig. 4), we have (e.g. Vidotto et al. 2014)

$$ \frac{r_M}{R_p} \approx f \left[ \frac{(B_p/2)^2/(8\pi)}{p_{\text{tot}}(a_{\text{orb}})} \right]^{1/6}, \quad (5) $$

where $B_p$ is the polar field strength of the planet's dipolar magnetic field (i.e. twice the value of the equatorial field strength) and $f \approx 2^{2/6}$ is a correction factor used to account for the effects of currents (e.g. Cravens 2004). The value of $B_p$ for GJ 436b is currently poorly constrained: Bourrier et al. (2016) derived an upper limit for the equatorial field strength of ~3 G, consistent with the upper limit estimates of ~10 G from Loyd et al. (2023). Assuming $B_p = 2$ G leads to an average magnetosphere size of 5.2 $R_p$ and 3.7 $R_p$ for stellar wind Models I and II, respectively, as seen on the top panel of Fig. 5. The maximum power that is radiated by Alfvén waves travelling towards the star is shown on the middle panel of this figure, where we see on average that the radiated power is $\mathcal{P}_{\text{wg}} \approx 1.2 \times 10^{22}$ erg s$^{-1}$ and $4.8 \times 10^{22}$ erg s$^{-1}$ for Models I and II, respectively. It is worth noting that the radiated power depends on the assumed magnetic field strength of the planet indirectly through $r_M$. As $r_M \propto B_p^{1/3}$, we have that $\mathcal{P}_{\text{wg}} \propto B_p^{2/3}$. Thus, for a 20 G polar field strength (i.e. 10 G at the equator), the radiated power increases by $10^{2/3}$, which is nearly a factor of 5, resulting in $\mathcal{P}_{\text{wg}} \approx 5.6 \times 10^{22}$ erg s$^{-1}$ and $22 \times 10^{22}$ erg s$^{-1}$ for Models I and II, respectively.

## 4.2. Scenario 2: Star-planet magnetic reconnection

Magnetic reconnection occurs when magnetic field lines of different polarities interact with each other. Magnetic reconnection has been considered in the case of star-planet interactions by several authors (e.g. Jardine & Cameron 2008; Lanza 2009; Vidotto et al. 2010). While Jardine & Cameron (2008); Vidotto et al. (2010) considered reconnection events to power radio emission from exoplanets, Lanza (2009) studied how reconnection events could trigger chromospheric hot spots. We can estimate the power released in the reconnection by integrating the Poynting flux, $P$ (Eq. (3)) over the cross-section of the planet's magnetosphere (Lanza 2009):

$$ \mathcal{P}_{\text{rec}} = \gamma_{\text{rec}} P(\pi R_M^2) = \frac{1}{4} \gamma_{\text{rec}} r_M^2 B^2 \Delta u \sin\theta, \quad (6) $$

where $\gamma_{\text{rec}}$ is an efficiency of the reconnection related to the angle between the stellar and planetary field lines. Even though in the GJ 436 system the planet will move through different orientations of the stellar wind magnetic field (e.g. going through

<figure><img src="image_6.png" /><figcaption>Fig. 5. Strength of the stellar wind-planet interactions computed in this work. Top: magnetospheric size of the planet, assuming a planetary magnetic field strength of 2 G at the pole (i.e., 1 G at the equator). Middle: maximum dissipated power of Alfvén waves travelling towards the star (Saur et al. 2013). Bottom: maximum dissipated power from the star-planet magnetic reconnection scenario (Lanza 2009). The solid lines are for Model I, and the dashed lines for Model II. The colour scheme follows the one adopted in Fig. 3. The grey arrows in the upper horizontal axes indicate approximate times of mid-transits.</figcaption></figure>

closed and open field lines), and thus $\gamma_{\text{rec}}$ would differ at different points in the orbit, for simplicity here we take its maximum value of $\gamma_{\text{rec}} = 1$, so that our power estimates $\mathcal{P}_{\text{rec}}$ are upper limits. The results of this scenario are shown in the bottom panel of Fig. 5, where we see that the average powers released in the reconnection scenario are $\mathcal{P}_{\text{rec}} \approx 8 \times 10^{22}$ erg s$^{-1}$ and $17 \times 10^{22}$ erg s$^{-1}$ for Models I and II, respectively, where we
assumed a planetary polar field strength of 2 G. Similarly to the previous scenario we considered, the power released in the reconnection scenario also depends on the planetary magnetic field as $P_{\text{rec}} \propto B_p^{2/3}$. Thus, for a 20 G field we find an increase in dissipated power by nearly a factor of 5; it reaches average values of $P_{\text{rec}} \approx 38 \times 10^{22}$ erg s$^{-1}$ and $80 \times 10^{22}$ erg s$^{-1}$ for Models I and II, respectively.

### 4.3. Interconnecting magnetic loop between star and planet

We calculated the power originating from star-planet interactions mediated by the stellar wind using the two scenarios discussed before. There is still a third model that we explore in the appendix, namely that of an interconnecting loop that extends from the stellar surface to the planetary orbit (Lanza 2013). For similar planetary field strengths, the interconnecting loop scenario produces maximum powers that are about three to four orders of magnitude larger than the other scenarios. This scenario assumes that the stellar magnetic field remains closed out to the orbital distance of the planet, which is not seen in our wind models.

For the three scenarios explored here (in Sect. 4 and in the appendix), the predicted powers show peaks and troughs, which are due to the planet crossing the magnetic equator of the star (when the field changes polarity). To compute the visibility of this interaction, one would need to take into account inclination of the rotation axis of the star (35.7$^{+5.9}_{-7.6}$ deg; Bourrier et al. 2022), the timing when the planet is hidden behind the stellar disk, and possible phase lags.

## 5. Discussion and conclusions

In this study we performed 3D numerical simulations of the wind of GJ 436. Our model incorporates the recently reconstructed radial magnetic field of GJ 436 (Bellotti et al. 2023), which shows that, in March-June 2016, GJ 436 presented a large-scale field that largely resembled an aligned dipole, with a maximum (absolute) field strength of ~27 G (Fig. 1). Stellar outflows, in the form of winds and coronal mass ejections, shape the space weather environment around planets. From our simulations, we derived the space weather conditions around GJ 436b, a warm-Neptune orbiting at a distance of 0.028 au from its host star. The main unknown in our wind models is the energy flux of the Alfvén waves that heat and accelerate the stellar wind. We ran two models with different values of the Alfvén wave fluxes to magnetic field ratio: the canonical value of $S_{A,0}/B_\star = 1.1 \times 10^7$ erg (cm$^2$ s G)$^{-1}$, commonly used in solar wind simulations, and another value that is one order of magnitude smaller. These models were named Model I (with the lower flux) and Model II (with the higher flux), as can be seen in Fig. 2. With the exception of the Alfvén wave properties, the remaining inputs of our models are well constrained by observations, including the large-scale surface field of the star. With our models, we derive stellar wind mass-loss rates of $1.1 \times 10^{-15}$ $M_\odot$ yr$^{-1}$ and $2.5 \times 10^{-14}$ $M_\odot$ yr$^{-1}$ for Models I and II, respectively. As discussed further below, we find that Model I is a better representation of the wind of GJ 436.

In addition to the global properties of the stellar wind, our models also allowed us to derive the conditions of the stellar wind at the orbit of GJ 436b. Even though at the observed epoch the stellar magnetic field, and thus the stellar wind, was approximately axi-symmetric about the rotation axis of the star, GJ 436b has a nearly polar orbit. This means that during one orbit, the

planet probes regions of open field lines and closed field lines. Taking the spin-orbit angle into account, we calculated the position of the planet during consecutive orbits: during one stellar rotation (~44 days) the planet orbits about ~16 times, experiencing different space weather conditions during each orbit (see Figs. 3 and 4). There are two main sources of variation in these conditions: one relates to the choice of wind model (lower versus higher Alfvén wave fluxes), and the other relates to the high obliquity of the orbit of GJ 436b. For example, the many peaks and troughs we see in the local wind properties (Fig. 4) are associated with the planet crossing the stellar equatorial plane.

Regardless of the two wind models we chose, our 3D simulations indicate that the motion of GJ 436b through the stellar wind is sub-Alfvénic, which implies that the planet has a direct magnetic connectivity with its host star. In a second part of our study, we investigated the power released in star-planet interactions mediated by the magnetised stellar wind. We first computed the Poynting fluxes of the stellar wind that impacts GJ 436b to be $P \approx 160$ and $650$ erg cm$^{-2}$ s$^{-1}$ (or 0.1 to 0.6 W m$^{-2}$) for Models I and II, respectively. Next, we computed the maximum energy that is available to power sub-Alfvénic star-planet interactions through two different scenarios: the Alfvén wing scenario (in which perturbations caused in the stellar magnetic field line by the planet travel towards the star as MHD waves, irradiating their energies; Saur et al. 2013) and the magnetic reconnection scenario (in which magnetic reconnection between stellar and planetary magnetic fields releases energy; Lanza 2009). For both scenarios, we assumed that the planet is magnetised with an equatorial dipolar magnetic field of 1G, leading to a magnetospheric size of about $5.2 R_p$ and $3.7 R_p$ for stellar wind Models I and II, respectively. Our estimated maximum powers released through star-planet interactions are of the order of $10^{22}$ to $10^{23}$ erg s$^{-1}$ (approximately $10^{-10}$ to $10^{-9}$ of the stellar bolometric luminosity), with the largest values associated with the reconnection scenario. It is interesting to note that the powers we derived are not very sensitive to the choice of Alfvén wave surface flux in the wind model (the largest unknown in our wind models) – a change of one order of magnitude in the surface flux of Alfvén waves leads to average powers that are a factor of 4 different for the Alfvén wing scenario and a factor of 2 in the reconnection scenario only (Fig. 5).

We note here that the planetary magnetic field is an assumption in our calculations and that the power released in the two star-planet interaction scenarios we considered increases with $B_p^{2/3}$. Thus, increasing the magnetic field of the planet by a factor of 10 (i.e. to an equatorial field of 10 G) leads to an increase in the estimated powers by $10^{2/3}$ (i.e. nearly a factor of 5). We could keep increasing the choice of planetary field strength until our estimated powers increased from ~$10^{23}$ erg s$^{-1}$ to $10^{26}$ erg s$^{-1}$, the latter of which is the maximum power derived by Loyd et al. (2023) for star-planet interactions in the GJ 436 system. To match the observed upper limit, the maximum planetary magnetic field would need to be unreasonably high ($\le 10^{3 \times 3/2} \approx 3 \times 10^4$ G). Because of this, we unfortunately cannot derive any meaningful upper limit to the planetary magnetic field.

According to the models studied here, GJ 436b orbits within the Alfvén surface (i.e. it has a sub-Alfvénic orbital motion) except for small portions of the orbit around the magnetic equator, where the orbit is outside the Alfvén surface. When a planet moves from sub- to super-Alfvénic motion, different types of star-planet interactions can take place. The formation of a bow shock around the planet should be seen as soon as the planet
assumed a planetary polar field strength of 2 G. Similarly to the previous scenario we considered, the power released in the reconnection scenario also depends on the planetary magnetic field as $P_{\text{rec}} \propto B_p^{2/3}$. Thus, for a 20 G field we find an increase in dissipated power by nearly a factor of 5; it reaches average values of $P_{\text{rec}} \approx 38 \times 10^{22}$ erg s$^{-1}$ and $80 \times 10^{22}$ erg s$^{-1}$ for Models I and II, respectively.

### 4.3. Interconnecting magnetic loop between star and planet

We calculated the power originating from star-planet interactions mediated by the stellar wind using the two scenarios discussed before. There is still a third model that we explore in the appendix, namely that of an interconnecting loop that extends from the stellar surface to the planetary orbit (Lanza 2013). For similar planetary field strengths, the interconnecting loop scenario produces maximum powers that are about three to four orders of magnitude larger than the other scenarios. This scenario assumes that the stellar magnetic field remains closed out to the orbital distance of the planet, which is not seen in our wind models.

For the three scenarios explored here (in Sect. 4 and in the appendix), the predicted powers show peaks and troughs, which are due to the planet crossing the magnetic equator of the star (when the field changes polarity). To compute the visibility of this interaction, one would need to take into account inclination of the rotation axis of the star (35.7$^{+5.9}_{-7.6}$ deg; Bourrier et al. 2022), the timing when the planet is hidden behind the stellar disk, and possible phase lags.

## 5. Discussion and conclusions

In this study we performed 3D numerical simulations of the wind of GJ 436. Our model incorporates the recently reconstructed radial magnetic field of GJ 436 (Bellotti et al. 2023), which shows that, in March-June 2016, GJ 436 presented a large-scale field that largely resembled an aligned dipole, with a maximum (absolute) field strength of ~27 G (Fig. 1). Stellar outflows, in the form of winds and coronal mass ejections, shape the space weather environment around planets. From our simulations, we derived the space weather conditions around GJ 436b, a warm-Neptune orbiting at a distance of 0.028 au from its host star. The main unknown in our wind models is the energy flux of the Alfvén waves that heat and accelerate the stellar wind. We ran two models with different values of the Alfvén wave fluxes to magnetic field ratio: the canonical value of $S_{A,0}/B_\star = 1.1 \times 10^7$ erg (cm$^2$ s G)$^{-1}$, commonly used in solar wind simulations, and another value that is one order of magnitude smaller. These models were named Model I (with the lower flux) and Model II (with the higher flux), as can be seen in Fig. 2. With the exception of the Alfvén wave properties, the remaining inputs of our models are well constrained by observations, including the large-scale surface field of the star. With our models, we derive stellar wind mass-loss rates of $1.1 \times 10^{-15}$ $M_\odot$ yr$^{-1}$ and $2.5 \times 10^{-14}$ $M_\odot$ yr$^{-1}$ for Models I and II, respectively. As discussed further below, we find that Model I is a better representation of the wind of GJ 436.

In addition to the global properties of the stellar wind, our models also allowed us to derive the conditions of the stellar wind at the orbit of GJ 436b. Even though at the observed epoch the stellar magnetic field, and thus the stellar wind, was approximately axi-symmetric about the rotation axis of the star, GJ 436b has a nearly polar orbit. This means that during one orbit, the

planet probes regions of open field lines and closed field lines. Taking the spin-orbit angle into account, we calculated the position of the planet during consecutive orbits: during one stellar rotation (~44 days) the planet orbits about ~16 times, experiencing different space weather conditions during each orbit (see Figs. 3 and 4). There are two main sources of variation in these conditions: one relates to the choice of wind model (lower versus higher Alfvén wave fluxes), and the other relates to the high obliquity of the orbit of GJ 436b. For example, the many peaks and troughs we see in the local wind properties (Fig. 4) are associated with the planet crossing the stellar equatorial plane.

Regardless of the two wind models we chose, our 3D simulations indicate that the motion of GJ 436b through the stellar wind is sub-Alfvénic, which implies that the planet has a direct magnetic connectivity with its host star. In a second part of our study, we investigated the power released in star-planet interactions mediated by the magnetised stellar wind. We first computed the Poynting fluxes of the stellar wind that impacts GJ 436b to be $P \approx 160$ and $650$ erg cm$^{-2}$ s$^{-1}$ (or 0.1 to 0.6 W m$^{-2}$) for Models I and II, respectively. Next, we computed the maximum energy that is available to power sub-Alfvénic star-planet interactions through two different scenarios: the Alfvén wing scenario (in which perturbations caused in the stellar magnetic field line by the planet travel towards the star as MHD waves, irradiating their energies; Saur et al. 2013) and the magnetic reconnection scenario (in which magnetic reconnection between stellar and planetary magnetic fields releases energy; Lanza 2009). For both scenarios, we assumed that the planet is magnetised with an equatorial dipolar magnetic field of 1G, leading to a magnetospheric size of about $5.2 R_p$ and $3.7 R_p$ for stellar wind Models I and II, respectively. Our estimated maximum powers released through star-planet interactions are of the order of $10^{22}$ to $10^{23}$ erg s$^{-1}$ (approximately $10^{-10}$ to $10^{-9}$ of the stellar bolometric luminosity), with the largest values associated with the reconnection scenario. It is interesting to note that the powers we derived are not very sensitive to the choice of Alfvén wave surface flux in the wind model (the largest unknown in our wind models) – a change of one order of magnitude in the surface flux of Alfvén waves leads to average powers that are a factor of 4 different for the Alfvén wing scenario and a factor of 2 in the reconnection scenario only (Fig. 5).

We note here that the planetary magnetic field is an assumption in our calculations and that the power released in the two star-planet interaction scenarios we considered increases with $B_p^{2/3}$. Thus, increasing the magnetic field of the planet by a factor of 10 (i.e. to an equatorial field of 10 G) leads to an increase in the estimated powers by $10^{2/3}$ (i.e. nearly a factor of 5). We could keep increasing the choice of planetary field strength until our estimated powers increased from ~$10^{23}$ erg s$^{-1}$ to $10^{26}$ erg s$^{-1}$, the latter of which is the maximum power derived by Loyd et al. (2023) for star-planet interactions in the GJ 436 system. To match the observed upper limit, the maximum planetary magnetic field would need to be unreasonably high ($\le 10^{3 \times 3/2} \approx 3 \times 10^4$ G). Because of this, we unfortunately cannot derive any meaningful upper limit to the planetary magnetic field.

According to the models studied here, GJ 436b orbits within the Alfvén surface (i.e. it has a sub-Alfvénic orbital motion) except for small portions of the orbit around the magnetic equator, where the orbit is outside the Alfvén surface. When a planet moves from sub- to super-Alfvénic motion, different types of star-planet interactions can take place. The formation of a bow shock around the planet should be seen as soon as the planet
moves from sub- to super-Alfvénic motion (e.g. Vidotto et al. 2011a). Cohen et al. (2014) studied the sub- to super-Alfvénic transition in terms of joule heating in the upper atmosphere of planets, while Kislyakova et al. (2018) and Kislyakova & Noack (2020) demonstrated that planets orbiting in the varying magnetic fields of their host stars can also experience induction heating in their interiors. Additionally, it is believed that magnetic connectivity can switch between on and off states – this can be due to variations over one planetary year (e.g. the planet moving in and out of the Alfvén surface) or variations in the Alfvén surface due to the evolution of the surface magnetic field (Fares et al. 2010). The latter fact has been attributed as the cause of the disappearance of planet-induced chromospheric modulation in planet-hosting stars (Shkolnik et al. 2008), emphasising the importance of contemporaneous observations of stellar magnetism and the search for signatures of star-planet interactions.

We explore a third star-planet interaction scenario in the appendix of this paper – in this scenario, a magnetic loop connecting the star and the planet would suffer a rupture after being stretched by the planet (Lanza 2013). Of the three scenarios, this has been recognised as the one that releases the largest star-planet interaction powers (Cauley et al. 2019; Strugarek et al. 2022). However, this scenario is not applicable to the case where the star-planet interaction is mediated by a stellar wind that distorts the stellar magnetic field from a potential state. In this case, in the calculations shown in the appendix, a potential field model for the stellar magnetic field is assumed (without a stellar wind). Using this scenario, we demonstrate in the appendix that a planetary field strength of ~6 G would reproduce the powers found by Loyd et al. (2023), supporting with these authors’ estimates.

In our work we adopted two different values of the Alfvén wave flux that drives the stellar wind. One way to constrain the best choice of surface Alfvén wave flux, and thus decide which of the wind models is more realistic, is to compare the mass-loss rate and wind properties we predict with the values obtained from models of the transmission spectroscopic transit of GJ 436b from Bourrier et al. (2016). Comparing the wind velocities and densities that we derived to those derived in Bourrier et al. (2016, shown as grey crosses in the upper panels of Fig. 4), we find that our Model I, with the lowest surface flux of Alfvén waves, matches the spectroscopic transit values better. Our derived mass-loss rate of $1.1 \times 10^{-15} M_{\odot} \text{ yr}^{-1}$ is compatible with that derived from spectroscopic transits of $1.2^{+1.3}_{-0.75} \times 10^{-15} M_{\odot} \text{ yr}^{-1}$ (Vidotto & Bourrier 2017). In addition to this, our preference towards Model I is also supported by further empirical and theoretical works, as we explain below.

Empirically, for the Sun, it has been suggested that the non-thermal velocities, $\xi$, of lines formed at the transition region are associated with the rms amplitude of the Alfvén waves (e.g. Oran et al. 2017) as $\xi = \langle v_{\perp 0}^2 \rangle^{1/2} / \sqrt{2}$. If this can be extended to solar-like stars, it implies that by measuring $\xi$ one can constrain the wave properties of other stars. Recently, Boro Saikia et al. (2023) measured the non-thermal velocity amplitudes of a sample of 55 low-mass stars. Overall, they find an increase in $\xi$ with stellar rotation (or activity). In their sample, they have five stars with masses, radii, and effective temperatures similar to those of GJ 436, namely: GJ 832, LTT 2050, GJ 3470, GJ 176, and GJ 849. Three of these stars (GJ 832, GJ 176, and GJ 849) have a rotation period of about 40 days, very similar to GJ 436. Using the Si IV line (1393 Å), Boro Saikia et al. (2023) measured $\xi$ ranging between 14.7 and 17.8 km s$^{-1}$ for these objects. These values are very similar to the non-thermal velocity of 14 km s$^{-1}$ from Model I.

Theoretically, it is usually assumed that the wave flux is given by a combination of the available turbulent energy flux inside a magnetic flux tube times the magnetic filling factor (e.g. Cranmer 2017; Suzuki 2018). The magnetic filling factor is smaller for slowly rotating stars (See et al. 2019), which are also usually older. Compared to the Sun, the turbulent energy flux is also lower in M dwarfs, given that the energy flux emerging from convective motions (and thus eventually driving the waves) is proportional to the bolometric energy flux of the star, $\sigma T_{\text{eff}}^4$ (e.g. Suzuki 2018). Sakaue & Shibata (2021) performed a parametric study of wave-driven wind models of M dwarfs, showing that the wave amplitude in the lower corona is smaller for M dwarfs than for the Sun.

In conclusion, both models (Sakaue & Shibata 2021) and empirical measurements (Boro Saikia et al. 2023) of the wave amplitude point to lower wave flux values for slowly rotating M dwarfs. Based on these arguments and the agreement of the local wind values we derive here with results from Ly-$\alpha$ transit modelling (Bourrier et al. 2016; Vidotto & Bourrier 2017), we favour a lower flux of Alfvén waves for the case of GJ 436. This might imply that the canonical value of $S_{A,0}/B_{\star} = 1.1 \times 10^7 \text{ erg (cm}^2 \text{ s G)}^{-1}$, derived from solar wind models, might not be applicable to the winds of old M dwarf stars.

Our models indicate that the planet GJ 436b experiences a vast range of stellar wind speeds and densities along its orbit. Because the wind shapes the atmospheric material that is being evaporated from the planet, we expect variations in the interaction region between the escaping atmosphere and the stellar wind, and hence, variations in the large structures formed due to this interaction, such as comet-like tails. We see the effects of different stellar wind properties on the morphology of comet-like tails, for example in 3D models that investigate these interactions (e.g. Villarreal D'Angelo et al. 2014, 2021; McCann et al. 2019; Carolan et al. 2021; Rumenskikh et al. 2022). Time-dependent effects, such as coronal mass ejections, would also alter the morphology of comet-like tails (Cherenkov et al. 2017; Hazra et al. 2022), similar to what we see when the solar wind interacts with the plasma tail of comets (Vourlidas et al. 2007). Because of the large range of stellar wind properties experienced by GJ 436b, similar temporal effects are expected during its orbit. It is interesting to note, however, the remarkable stability of the wind properties at every consecutive mid-transit time (crosses in the upper panel of Fig. 4). If the local wind properties are relatively constant for a few hours, this would lead to a relatively stable interaction zone, potentially explaining why the Ly-$\alpha$ transit of GJ 436b remains stable over time (Dos Santos et al. 2019).

For this stability to occur in transits observed many years apart, the stellar magnetic field cannot vary substantially. To investigate this, we need further spectropolarimetric campaigns to reconstruct the large-scale field of the star and thus study its evolution. We note that GJ 436 has an activity cycle with a period of about 7–8 yr (Lothringer et al. 2018; Dos Santos et al. 2019; Kumar & Fares 2023; Loyd et al. 2023), which implies that the stellar magnetic field could vary (in strength and topology) substantially at different phases of the cycle. This temporal variation would affect the wind properties and space weather conditions around the planet (e.g. Vidotto et al. 2012; Nicholson et al. 2016; Finley et al. 2018; Kavanagh et al. 2019) and, consequently, the power released during star-planet interactions (e.g. Klein et al. 2022). Ideally, observations that probe star-planet interactions mediated by the stellar wind and/or magnetic field (e.g. planetary radio emission, chromospheric hot spots, or the interaction of
moves from sub- to super-Alfvénic motion (e.g. Vidotto et al. 2011a). Cohen et al. (2014) studied the sub- to super-Alfvénic transition in terms of joule heating in the upper atmosphere of planets, while Kislyakova et al. (2018) and Kislyakova & Noack (2020) demonstrated that planets orbiting in the varying magnetic fields of their host stars can also experience induction heating in their interiors. Additionally, it is believed that magnetic connectivity can switch between on and off states – this can be due to variations over one planetary year (e.g. the planet moving in and out of the Alfvén surface) or variations in the Alfvén surface due to the evolution of the surface magnetic field (Fares et al. 2010). The latter fact has been attributed as the cause of the disappearance of planet-induced chromospheric modulation in planet-hosting stars (Shkolnik et al. 2008), emphasising the importance of contemporaneous observations of stellar magnetism and the search for signatures of star-planet interactions.

We explore a third star-planet interaction scenario in the appendix of this paper – in this scenario, a magnetic loop connecting the star and the planet would suffer a rupture after being stretched by the planet (Lanza 2013). Of the three scenarios, this has been recognised as the one that releases the largest star-planet interaction powers (Cauley et al. 2019; Strugarek et al. 2022). However, this scenario is not applicable to the case where the star-planet interaction is mediated by a stellar wind that distorts the stellar magnetic field from a potential state. In this case, in the calculations shown in the appendix, a potential field model for the stellar magnetic field is assumed (without a stellar wind). Using this scenario, we demonstrate in the appendix that a planetary field strength of ~6 G would reproduce the powers found by Loyd et al. (2023), supporting with these authors’ estimates.

In our work we adopted two different values of the Alfvén wave flux that drives the stellar wind. One way to constrain the best choice of surface Alfvén wave flux, and thus decide which of the wind models is more realistic, is to compare the mass-loss rate and wind properties we predict with the values obtained from models of the transmission spectroscopic transit of GJ 436b from Bourrier et al. (2016). Comparing the wind velocities and densities that we derived to those derived in Bourrier et al. (2016, shown as grey crosses in the upper panels of Fig. 4), we find that our Model I, with the lowest surface flux of Alfvén waves, matches the spectroscopic transit values better. Our derived mass-loss rate of $1.1 \times 10^{-15} M_{\odot} \text{ yr}^{-1}$ is compatible with that derived from spectroscopic transits of $1.2^{+1.3}_{-0.75} \times 10^{-15} M_{\odot} \text{ yr}^{-1}$ (Vidotto & Bourrier 2017). In addition to this, our preference towards Model I is also supported by further empirical and theoretical works, as we explain below.

Empirically, for the Sun, it has been suggested that the non-thermal velocities, $\xi$, of lines formed at the transition region are associated with the rms amplitude of the Alfvén waves (e.g. Oran et al. 2017) as $\xi = \langle v_{\perp 0}^2 \rangle^{1/2} / \sqrt{2}$. If this can be extended to solar-like stars, it implies that by measuring $\xi$ one can constrain the wave properties of other stars. Recently, Boro Saikia et al. (2023) measured the non-thermal velocity amplitudes of a sample of 55 low-mass stars. Overall, they find an increase in $\xi$ with stellar rotation (or activity). In their sample, they have five stars with masses, radii, and effective temperatures similar to those of GJ 436, namely: GJ 832, LTT 2050, GJ 3470, GJ 176, and GJ 849. Three of these stars (GJ 832, GJ 176, and GJ 849) have a rotation period of about 40 days, very similar to GJ 436. Using the Si IV line (1393 Å), Boro Saikia et al. (2023) measured $\xi$ ranging between 14.7 and 17.8 km s$^{-1}$ for these objects. These values are very similar to the non-thermal velocity of 14 km s$^{-1}$ from Model I.

Theoretically, it is usually assumed that the wave flux is given by a combination of the available turbulent energy flux inside a magnetic flux tube times the magnetic filling factor (e.g. Cranmer 2017; Suzuki 2018). The magnetic filling factor is smaller for slowly rotating stars (See et al. 2019), which are also usually older. Compared to the Sun, the turbulent energy flux is also lower in M dwarfs, given that the energy flux emerging from convective motions (and thus eventually driving the waves) is proportional to the bolometric energy flux of the star, $\sigma T_{\text{eff}}^4$ (e.g. Suzuki 2018). Sakaue & Shibata (2021) performed a parametric study of wave-driven wind models of M dwarfs, showing that the wave amplitude in the lower corona is smaller for M dwarfs than for the Sun.

In conclusion, both models (Sakaue & Shibata 2021) and empirical measurements (Boro Saikia et al. 2023) of the wave amplitude point to lower wave flux values for slowly rotating M dwarfs. Based on these arguments and the agreement of the local wind values we derive here with results from Ly-$\alpha$ transit modelling (Bourrier et al. 2016; Vidotto & Bourrier 2017), we favour a lower flux of Alfvén waves for the case of GJ 436. This might imply that the canonical value of $S_{A,0}/B_{\star} = 1.1 \times 10^7 \text{ erg (cm}^2 \text{ s G)}^{-1}$, derived from solar wind models, might not be applicable to the winds of old M dwarf stars.

Our models indicate that the planet GJ 436b experiences a vast range of stellar wind speeds and densities along its orbit. Because the wind shapes the atmospheric material that is being evaporated from the planet, we expect variations in the interaction region between the escaping atmosphere and the stellar wind, and hence, variations in the large structures formed due to this interaction, such as comet-like tails. We see the effects of different stellar wind properties on the morphology of comet-like tails, for example in 3D models that investigate these interactions (e.g. Villarreal D'Angelo et al. 2014, 2021; McCann et al. 2019; Carolan et al. 2021; Rumenskikh et al. 2022). Time-dependent effects, such as coronal mass ejections, would also alter the morphology of comet-like tails (Cherenkov et al. 2017; Hazra et al. 2022), similar to what we see when the solar wind interacts with the plasma tail of comets (Vourlidas et al. 2007). Because of the large range of stellar wind properties experienced by GJ 436b, similar temporal effects are expected during its orbit. It is interesting to note, however, the remarkable stability of the wind properties at every consecutive mid-transit time (crosses in the upper panel of Fig. 4). If the local wind properties are relatively constant for a few hours, this would lead to a relatively stable interaction zone, potentially explaining why the Ly-$\alpha$ transit of GJ 436b remains stable over time (Dos Santos et al. 2019).

For this stability to occur in transits observed many years apart, the stellar magnetic field cannot vary substantially. To investigate this, we need further spectropolarimetric campaigns to reconstruct the large-scale field of the star and thus study its evolution. We note that GJ 436 has an activity cycle with a period of about 7–8 yr (Lothringer et al. 2018; Dos Santos et al. 2019; Kumar & Fares 2023; Loyd et al. 2023), which implies that the stellar magnetic field could vary (in strength and topology) substantially at different phases of the cycle. This temporal variation would affect the wind properties and space weather conditions around the planet (e.g. Vidotto et al. 2012; Nicholson et al. 2016; Finley et al. 2018; Kavanagh et al. 2019) and, consequently, the power released during star-planet interactions (e.g. Klein et al. 2022). Ideally, observations that probe star-planet interactions mediated by the stellar wind and/or magnetic field (e.g. planetary radio emission, chromospheric hot spots, or the interaction of
evaporated atmospheres with stellar winds) and the magnetic map used in stellar wind models should be conducted contemporaneously.

**Acknowledgements.** A.A.V. thanks P. Loyd and A. Lanza for useful discussions in the preparation of this paper. The authors thank the anonymous reviewer for their constructive recommendations and suggestions. This project has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (grant agreement No 817540, ASTROFLOW; grant agreement No 947634, SPICE DUNE; grant agreement No 740651, NewWorlds). RF acknowledges support from the United Arab Emirates University (UAEU) startup grant number G00003269. This work has been carried out within the framework of the NCCR PlanetS supported by the Swiss National Science Foundation under grants 51NF40_182901 and 51NF40_205606. We acknowledge funding from the French National Research Agency (ANR) under contract number ANR-18-CE31-0019 (SPlaSH). We thank SURF (www.surf.nl) for the support in using the National Supercomputer Snellius. This work used the BATS-R-US tools developed at the University of Michigan Center for Space Environment Modeling and made available through the NASA Community Coordinated Modeling Center.

## References

Adams, F. C. 2011, ApJ, 730, 27

Attia, O., Bourrier, V., Eggenberger, P., et al. 2021, A&A, 647, A40

Bellotti, S., Fares, R., Vidotto, A. A., et al. 2023, A&A, 676, A139

Boro Saikia, S., Jin, M., Johnstone, C. P., et al. 2020, A&A, 635, A178

Boro Saikia, S., Lueftinger, T., Airapetian, V. S., et al. 2023, ApJ, 950, 124

Bourrier, V., Lecavelier des Etangs, A., Ehrenreich, D., Tanaka, Y. A., & Vidotto, A. A. 2016, A&A, 591, A121

Bourrier, V., Lovis, C., Beust, H., et al. 2018, Nature, 553, 477

Bourrier, V., Zapatero Osorio, M. R., Allart, R., et al. 2022, A&A, 663, A160

Carolan, S., Vidotto, A. A., Villarreal D'Angelo, C., & Hazra, G. 2021, MNRAS, 500, 3382

Cauley, P. W., Shkolnik, E. L., Llama, J., & Lanza, A. F. 2019, Nat. Astron., 3, 1128

Cherenkov, A., Bisikalo, D., Fossati, L., & Mostl, C. 2017, ApJ, 846, 31

Cohen, O., Drake, J. J., Glocer, A., et al. 2014, ApJ, 790, 57

Cranmer, S. R. 2017, ApJ, 840, 114

Cravens, T. E. 2004, Physics of Solar System Plasmas (Cambridge: Cambridge University Press)

De Pontieu, B., McIntosh, S. W., Carlsson, M., et al. 2007, Science, 318, 1574

Donati, J., & Landstreet, J. D. 2009, ARA&A, 47, 333

Donati, J.-F., & Brown, S. F. 1997, A&A, 326, 1135

Dos Santos, L. A., Ehrenreich, D., Bourrier, V., et al. 2019, A&A, 629, A47

Ehrenreich, D., Bourrier, V., Wheatley, P. J., et al. 2015, Nature, 522, 459

Elekes, F., & Saur, J. 2023, A&A, 671, A133

Evensberget, D., Marsden, S. C., Carter, B. D., et al. 2023, MNRAS, 524, 2042

Fares, R., Donati, J., Moutou, C., et al. 2010, MNRAS, 406, 409

Farrell, W. M., Desch, M. D., & Zarka, P. 1999, J. Geophys. Res., 104, 14025

Finley, A. J., Matt, S. P., & See, V. 2018, ApJ, 864, 125

Fischer, C., & Saur, J. 2019, ApJ, 872, 113

Folsom, C. P., O'Fionnagain, D., Fossati, L., et al. 2020, A&A, 633, A48

Fontenla, J. M., Linsky, J. L., Witbrod, J., et al. 2016, ApJ, 830, 154

Hazra, G., Vidotto, A. A., Carolan, S., Villarreal D'Angelo, C., & Manchester, W. 2022, MNRAS, 509, 5858

Hollweg, J. V. 1986, J. Geophys. Res., 91, 4111

Jardine, M., & Cameron, A. C. 2008, A&A, 490, 843

Kavanagh, R. D., Vidotto, A. A., O'Fionnagain, D., et al. 2019, MNRAS, 485, 4529

Kavanagh, R. D., Vidotto, A. A., Klein, B., et al. 2021, MNRAS, 504, 1511

Kavanagh, R. D., Vidotto, A. A., Vedantham, H. K., et al. 2022, MNRAS, 514, 675

Khodachenko, M. L., Shaikhislamov, I. F., Lammer, H., et al. 2019, ApJ, 885, 67

Kislyakova, K., & Noack, L. 2020, A&A, 636, L10

Kislyakova, K. G., Fossati, L., Johnstone, C. P., et al. 2018, ApJ, 858, 105

Klein, B., Zicher, N., Kavanagh, R. D., et al. 2022, MNRAS, 512, 5067

Kulow, J. R., France, K., Linsky, J., & Loyd, R. O. P. 2014, ApJ, 786, 132

Kumar, M., & Fares, R. 2023, MNRAS, 518, 3147

Lanotte, A. A., Gillon, M., Demory, B. O., et al. 2014, A&A, 572, A73

Lanza, A. F. 2009, A&A, 505, 339

Lanza, A. F. 2013, A&A, 557, A31

Lavie, B., Ehrenreich, D., Bourrier, V., et al. 2017, A&A, 605, L7

Lothringer, J. D., Benneke, B., Crossfield, I. J. M., et al. 2018, AJ, 155, 66

Loyd, R. O. P., Schneider, P. C., Jackman, J. A. G., et al. 2023, AJ, 165, 146

Mauas, P. J. D., Falchi, A., Pasquini, L., & Pallavicini, R. 1997, A&A, 326, 249

Maxted, P. F. L., Ehrenreich, D., Wilson, T. G., et al. 2022, MNRAS, 514, 77

McCann, J., Murray-Clay, R. A., Kratter, K., & Krumholz, M. R. 2019, ApJ, 873, 89

Mesquita, A. L., & Vidotto, A. A. 2020, MNRAS, 494, 1297

Mesquita, A. L., Rodgers-Lee, D., & Vidotto, A. A. 2021, MNRAS, 505, 1817

Nicholson, B. A., Vidotto, A. A., Mengel, M., et al. 2016, MNRAS, 459, 1907

O'Fionnagain, D., Vidotto, A. A., Petit, P., et al. 2021, MNRAS, 500, 3438

Oran, R., van der Holst, B., Landi, E., et al. 2013, ApJ, 778, 176

Oran, R., Landi, E., van der Holst, B., Sokolov, I. V., & Gombosi, T. I. 2017, ApJ, 845, 98

Pagano, I., Lanza, A. F., Leto, G., et al. 2009, Earth Moon Planets, 105, 373

Peacock, S., Barman, T., Shkolnik, E. L., et al. 2019, ApJ, 886, 77

Pineda, J. S., & Villadsen, J. 2023, Nat. Astron., 7, 569

Powell, K. G., Roe, P. L., Linde, T. J., Gombosi, T. I., & de Zeeuw, D. L. 1999, J. Chem. Phys., 154, 284

Preusse, S., Kopp, A., Buechner, J., & Motschmann, U. 2006, A&A, 460, 317

Rodgers-Lee, D., Rimmer, P. B., Vidotto, A. A., et al. 2023, MNRAS, 521, 5880

Rumenskikh, M. S., Shaikhislamov, I. F., Khodachenko, M. L., et al. 2022, ApJ, 927, 238

Sakaue, T., & Shibata, K. 2021, ApJ, 919, 29

Saur, J., Grambusch, T., Duling, S., Neubauer, F. M., & Simon, S. 2013, A&A, 552, A119

See, V., Matt, S. P., Folsom, C. P., et al. 2019, ApJ, 876, 118

Shaikhislamov, I. F., Khodachenko, M. L., Lammer, H., et al. 2018, MNRAS, 481, 5315

Shkolnik, E., Walker, G. A. H., Bohlender, D. A., Gu, P.-G., & Kuerster, M. 2005, ApJ, 622, 1075

Shkolnik, E., Bohlender, D. A., Walker, G. A. H., & Collier Cameron, A. 2008, ApJ, 676, 628

Shoda, M., Yokoyama, T., & Suzuki, T. K. 2018, ApJ, 853, 190

Sokolov, I. V., van der Holst, B., Oran, R., et al. 2013, ApJ, 764, 23

Strugarek, A. 2016, ApJ, 833, 140

Strugarek, A., Fares, R., Bourrier, V., et al. 2022, MNRAS, 512, 4556

Suzuki, T. K. 2018, PASJ, 70, 34

Torres, G. 2007, ApJ, 671, L65

Toth, G., van der Holst, B., Sokolov, I. V., et al. 2012, J. Comput. Phys., 231, 870

Trifonov, T., Kürster, M., Zechmeister, M., et al. 2018, A&A, 609, A117

Turner, J. D., Zarka, P., Grießmeier, J.-M., et al. 2021, A&A, 645, A59

van der Holst, B., Sokolov, I. V., Meng, X., et al. 2014, ApJ, 782, 81

Vidotto, A. A. 2020, IAU Symp., 354, 259

Vidotto, A. A. 2021, Liv. Rev. Sol. Phys., 18, 3

Vidotto, A. A., & Bourrier, V. 2017, MNRAS, 470, 4026

Vidotto, A. A., Opher, M., Jatenco-Pereira, V., & Gombosi, T. I. 2010, ApJ, 720, 1262

Vidotto, A. A., Jardine, M., & Helling, C. 2011a, MNRAS, 414, 1573

Vidotto, A. A., Jardine, M., Opher, M., Donati, J. F., & Gombosi, T. I. 2011b, MNRAS, 412, 351

Vidotto, A. A., Fares, R., Jardine, M., et al. 2012, MNRAS, 423, 3285

Vidotto, A. A., Jardine, M., Morin, J., et al. 2014, MNRAS, 438, 1162

Vidotto, A. A., Feeney, N., & Groh, J. H. 2019, MNRAS, 488, 633

Villarreal D'Angelo, C., Schneiter, M., Costa, A., et al. 2014, MNRAS, 438, 1654

Villarreal D'Angelo, C., Vidotto, A. A., Esquivel, A., Hazra, G., & Youngblood, A. 2021, MNRAS, 501, 4383

Vourlidas, A., Davis, C. J., Eyles, C. J., et al. 2007, ApJ, 668, L79
evaporated atmospheres with stellar winds) and the magnetic map used in stellar wind models should be conducted contemporaneously.

**Acknowledgements.** A.A.V. thanks P. Loyd and A. Lanza for useful discussions in the preparation of this paper. The authors thank the anonymous reviewer for their constructive recommendations and suggestions. This project has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (grant agreement No 817540, ASTROFLOW; grant agreement No 947634, SPICE DUNE; grant agreement No 740651, NewWorlds). RF acknowledges support from the United Arab Emirates University (UAEU) startup grant number G00003269. This work has been carried out within the framework of the NCCR PlanetS supported by the Swiss National Science Foundation under grants 51NF40_182901 and 51NF40_205606. We acknowledge funding from the French National Research Agency (ANR) under contract number ANR-18-CE31-0019 (SPlaSH). We thank SURF (www.surf.nl) for the support in using the National Supercomputer Snellius. This work used the BATS-R-US tools developed at the University of Michigan Center for Space Environment Modeling and made available through the NASA Community Coordinated Modeling Center.

## References

Adams, F. C. 2011, ApJ, 730, 27

Attia, O., Bourrier, V., Eggenberger, P., et al. 2021, A&A, 647, A40

Bellotti, S., Fares, R., Vidotto, A. A., et al. 2023, A&A, 676, A139

Boro Saikia, S., Jin, M., Johnstone, C. P., et al. 2020, A&A, 635, A178

Boro Saikia, S., Lueftinger, T., Airapetian, V. S., et al. 2023, ApJ, 950, 124

Bourrier, V., Lecavelier des Etangs, A., Ehrenreich, D., Tanaka, Y. A., & Vidotto, A. A. 2016, A&A, 591, A121

Bourrier, V., Lovis, C., Beust, H., et al. 2018, Nature, 553, 477

Bourrier, V., Zapatero Osorio, M. R., Allart, R., et al. 2022, A&A, 663, A160

Carolan, S., Vidotto, A. A., Villarreal D'Angelo, C., & Hazra, G. 2021, MNRAS, 500, 3382

Cauley, P. W., Shkolnik, E. L., Llama, J., & Lanza, A. F. 2019, Nat. Astron., 3, 1128

Cherenkov, A., Bisikalo, D., Fossati, L., & Mostl, C. 2017, ApJ, 846, 31

Cohen, O., Drake, J. J., Glocer, A., et al. 2014, ApJ, 790, 57

Cranmer, S. R. 2017, ApJ, 840, 114

Cravens, T. E. 2004, Physics of Solar System Plasmas (Cambridge: Cambridge University Press)

De Pontieu, B., McIntosh, S. W., Carlsson, M., et al. 2007, Science, 318, 1574

Donati, J., & Landstreet, J. D. 2009, ARA&A, 47, 333

Donati, J.-F., & Brown, S. F. 1997, A&A, 326, 1135

Dos Santos, L. A., Ehrenreich, D., Bourrier, V., et al. 2019, A&A, 629, A47

Ehrenreich, D., Bourrier, V., Wheatley, P. J., et al. 2015, Nature, 522, 459

Elekes, F., & Saur, J. 2023, A&A, 671, A133

Evensberget, D., Marsden, S. C., Carter, B. D., et al. 2023, MNRAS, 524, 2042

Fares, R., Donati, J., Moutou, C., et al. 2010, MNRAS, 406, 409

Farrell, W. M., Desch, M. D., & Zarka, P. 1999, J. Geophys. Res., 104, 14025

Finley, A. J., Matt, S. P., & See, V. 2018, ApJ, 864, 125

Fischer, C., & Saur, J. 2019, ApJ, 872, 113

Folsom, C. P., O'Fionnagain, D., Fossati, L., et al. 2020, A&A, 633, A48

Fontenla, J. M., Linsky, J. L., Witbrod, J., et al. 2016, ApJ, 830, 154

Hazra, G., Vidotto, A. A., Carolan, S., Villarreal D'Angelo, C., & Manchester, W. 2022, MNRAS, 509, 5858

Hollweg, J. V. 1986, J. Geophys. Res., 91, 4111

Jardine, M., & Cameron, A. C. 2008, A&A, 490, 843

Kavanagh, R. D., Vidotto, A. A., O'Fionnagain, D., et al. 2019, MNRAS, 485, 4529

Kavanagh, R. D., Vidotto, A. A., Klein, B., et al. 2021, MNRAS, 504, 1511

Kavanagh, R. D., Vidotto, A. A., Vedantham, H. K., et al. 2022, MNRAS, 514, 675

Khodachenko, M. L., Shaikhislamov, I. F., Lammer, H., et al. 2019, ApJ, 885, 67

Kislyakova, K., & Noack, L. 2020, A&A, 636, L10

Kislyakova, K. G., Fossati, L., Johnstone, C. P., et al. 2018, ApJ, 858, 105

Klein, B., Zicher, N., Kavanagh, R. D., et al. 2022, MNRAS, 512, 5067

Kulow, J. R., France, K., Linsky, J., & Loyd, R. O. P. 2014, ApJ, 786, 132

Kumar, M., & Fares, R. 2023, MNRAS, 518, 3147

Lanotte, A. A., Gillon, M., Demory, B. O., et al. 2014, A&A, 572, A73

Lanza, A. F. 2009, A&A, 505, 339

Lanza, A. F. 2013, A&A, 557, A31

Lavie, B., Ehrenreich, D., Bourrier, V., et al. 2017, A&A, 605, L7

Lothringer, J. D., Benneke, B., Crossfield, I. J. M., et al. 2018, AJ, 155, 66

Loyd, R. O. P., Schneider, P. C., Jackman, J. A. G., et al. 2023, AJ, 165, 146

Mauas, P. J. D., Falchi, A., Pasquini, L., & Pallavicini, R. 1997, A&A, 326, 249

Maxted, P. F. L., Ehrenreich, D., Wilson, T. G., et al. 2022, MNRAS, 514, 77

McCann, J., Murray-Clay, R. A., Kratter, K., & Krumholz, M. R. 2019, ApJ, 873, 89

Mesquita, A. L., & Vidotto, A. A. 2020, MNRAS, 494, 1297

Mesquita, A. L., Rodgers-Lee, D., & Vidotto, A. A. 2021, MNRAS, 505, 1817

Nicholson, B. A., Vidotto, A. A., Mengel, M., et al. 2016, MNRAS, 459, 1907

O'Fionnagain, D., Vidotto, A. A., Petit, P., et al. 2021, MNRAS, 500, 3438

Oran, R., van der Holst, B., Landi, E., et al. 2013, ApJ, 778, 176

Oran, R., Landi, E., van der Holst, B., Sokolov, I. V., & Gombosi, T. I. 2017, ApJ, 845, 98

Pagano, I., Lanza, A. F., Leto, G., et al. 2009, Earth Moon Planets, 105, 373

Peacock, S., Barman, T., Shkolnik, E. L., et al. 2019, ApJ, 886, 77

Pineda, J. S., & Villadsen, J. 2023, Nat. Astron., 7, 569

Powell, K. G., Roe, P. L., Linde, T. J., Gombosi, T. I., & de Zeeuw, D. L. 1999, J. Chem. Phys., 154, 284

Preusse, S., Kopp, A., Buechner, J., & Motschmann, U. 2006, A&A, 460, 317

Rodgers-Lee, D., Rimmer, P. B., Vidotto, A. A., et al. 2023, MNRAS, 521, 5880

Rumenskikh, M. S., Shaikhislamov, I. F., Khodachenko, M. L., et al. 2022, ApJ, 927, 238

Sakaue, T., & Shibata, K. 2021, ApJ, 919, 29

Saur, J., Grambusch, T., Duling, S., Neubauer, F. M., & Simon, S. 2013, A&A, 552, A119

See, V., Matt, S. P., Folsom, C. P., et al. 2019, ApJ, 876, 118

Shaikhislamov, I. F., Khodachenko, M. L., Lammer, H., et al. 2018, MNRAS, 481, 5315

Shkolnik, E., Walker, G. A. H., Bohlender, D. A., Gu, P.-G., & Kuerster, M. 2005, ApJ, 622, 1075

Shkolnik, E., Bohlender, D. A., Walker, G. A. H., & Collier Cameron, A. 2008, ApJ, 676, 628

Shoda, M., Yokoyama, T., & Suzuki, T. K. 2018, ApJ, 853, 190

Sokolov, I. V., van der Holst, B., Oran, R., et al. 2013, ApJ, 764, 23

Strugarek, A. 2016, ApJ, 833, 140

Strugarek, A., Fares, R., Bourrier, V., et al. 2022, MNRAS, 512, 4556

Suzuki, T. K. 2018, PASJ, 70, 34

Torres, G. 2007, ApJ, 671, L65

Toth, G., van der Holst, B., Sokolov, I. V., et al. 2012, J. Comput. Phys., 231, 870

Trifonov, T., Kürster, M., Zechmeister, M., et al. 2018, A&A, 609, A117

Turner, J. D., Zarka, P., Grießmeier, J.-M., et al. 2021, A&A, 645, A59

van der Holst, B., Sokolov, I. V., Meng, X., et al. 2014, ApJ, 782, 81

Vidotto, A. A. 2020, IAU Symp., 354, 259

Vidotto, A. A. 2021, Liv. Rev. Sol. Phys., 18, 3

Vidotto, A. A., & Bourrier, V. 2017, MNRAS, 470, 4026

Vidotto, A. A., Opher, M., Jatenco-Pereira, V., & Gombosi, T. I. 2010, ApJ, 720, 1262

Vidotto, A. A., Jardine, M., & Helling, C. 2011a, MNRAS, 414, 1573

Vidotto, A. A., Jardine, M., Opher, M., Donati, J. F., & Gombosi, T. I. 2011b, MNRAS, 412, 351

Vidotto, A. A., Fares, R., Jardine, M., et al. 2012, MNRAS, 423, 3285

Vidotto, A. A., Jardine, M., Morin, J., et al. 2014, MNRAS, 438, 1162

Vidotto, A. A., Feeney, N., & Groh, J. H. 2019, MNRAS, 488, 633

Villarreal D'Angelo, C., Schneiter, M., Costa, A., et al. 2014, MNRAS, 438, 1654

Villarreal D'Angelo, C., Vidotto, A. A., Esquivel, A., Hazra, G., & Youngblood, A. 2021, MNRAS, 501, 4383

Vourlidas, A., Davis, C. J., Eyles, C. J., et al. 2007, ApJ, 668, L79
## Appendix A: Calculation of orbital trajectories for a misaligned orbit

To calculate the orbital trajectories in the reference frame of the star, we compute two coordinate transformations: one from the orbital reference frame to the observer's frame, and another from the observer's frame to the stellar frame. These frames are illustrated in Fig. A.1.

We start from the orbital reference frame {$x_{pl}, y_{pl}, z_{pl}$}, such that the $x_{pl}y_{pl}$-plane contains the orbit of the planet and the orbital axis is along $z_{pl}$. The origin of this reference frame is centred on the star. In this reference frame, the (circular) planetary motion is simply described as

$$ \begin{bmatrix} x_{pl} \\ y_{pl} \\ z_{pl} \end{bmatrix} = \begin{bmatrix} R_{orb} \cos \varphi_{pl} \\ R_{orb} \sin \varphi_{pl} \\ 0 \end{bmatrix}, \qquad (A.1) $$

where the phase of the orbit is $\varphi_{pl} = \Omega_{orb}t + \varphi_0$, with $t$ describing the time and $\Omega_{orb} = 2\pi/P_{orb}$ the orbital rotation rate. Here, we assume a non-null initial phase $\varphi_0$, which we computed later in this section to find the location of GJ 436b at zero phase of the surface magnetic map, which corresponds to the first spectropolarimetric observation as reported in Bellotti et al. (2023).

Now we take the observer's (inertial) reference frame {$x_{obs}, y_{obs}, z_{obs}$}, whose origin is also centred on the star, and $z_{obs}$ is aligned with the stellar rotation. The observer is located at positive $x_{obs}$ and we assume that $x_{obs} \parallel x_{pl}$. The angle between the orbital spin and the stellar spin axis is the spin-orbit angle, $\Psi$; thus, the angle between $z_{obs}$ and $z_{pl}$ is also $\Psi$. In this reference frame, the orbital motion of the planet is described as

$$ \begin{bmatrix} x_{obs} \\ y_{obs} \\ z_{obs} \end{bmatrix} = \begin{bmatrix} 1 & 0 & 0 \\ 0 & \cos \Psi & \sin \Psi \\ 0 & -\sin \Psi & \cos \Psi \end{bmatrix} \begin{bmatrix} x_{pl} \\ y_{pl} \\ z_{pl} \end{bmatrix}, \qquad (A.2) $$

where we performed an anti-clockwise rotation around the axis $x_{obs} \equiv x_{pl}$ by an angle $\Psi$. According to our reference frame, the mid-transit occurs when the planet is at {$x_{obs}, y_{obs}, z_{obs}$} = {$R_{orb}, 0, 0$}.

<figure><img src="image_6.png" /><figcaption>Fig. A.1. Reference frames used in the coordinate transformation, with the star given by the orange circle, the planet indicated by the red circle, and its orbital path shown in red. The planet's reference frame {$x_{pl}, y_{pl}, z_{pl}$} (in red) contains the orbital plane of the planet $x_{pl}y_{pl}$ and orbital spin axis along $z_{pl}$. The observer's inertial reference frame {$x_{obs}, y_{obs}, z_{obs}$} (in blue) is centred on the star, with the observer located at $+x_{obs}$ and the stellar rotation axis aligned along $z_{obs}$. The reference frame {$x_{\star}, y_{\star}, z_{\star}$} (in black) co-rotates with the star. The spin-orbit angle is $\Psi$.</figcaption></figure>

The reference frame co-rotating with the star {$x_{\star}, y_{\star}, z_{\star}$} is also centred on the star, with the stellar spin axis along $z_{\star}$, which is parallel to $z_{obs}$. To transform from the observer's reference frame to the stellar co-rotating frame, we performed a clockwise rotation around the axis $z_{\star} \equiv z_{obs}$ by an angle $\Omega_{\star}t$:

$$ \begin{bmatrix} x_{\star} \\ y_{\star} \\ z_{\star} \end{bmatrix} = \begin{bmatrix} \cos(\Omega_{\star}t) & -\sin(\Omega_{\star}t) & 0 \\ \sin(\Omega_{\star}t) & \cos(\Omega_{\star}t) & 0 \\ 0 & 0 & 1 \end{bmatrix} \begin{bmatrix} x_{obs} \\ y_{obs} \\ z_{obs} \end{bmatrix}. \qquad (A.3) $$

Therefore, to derive the orbital path of the planet in the reference frame of the star, we substitute (A.1) and (A.2) into (A.3), thus obtaining

$$ \frac{x_{\star}}{R_{\text{orb}}} = \cos(\Omega_{\text{orb}}t + \varphi_0) \cos(\Omega_{\star}t) - \sin(\Omega_{\text{orb}}t + \varphi_0) \sin(\Omega_{\star}t) \cos \Psi \qquad (A.4) $$

$$ \frac{y_{\star}}{R_{\text{orb}}} = \cos(\Omega_{\text{orb}}t + \varphi_0) \sin(\Omega_{\star}t) + \sin(\Omega_{\text{orb}}t + \varphi_0) \cos(\Omega_{\star}t) \cos \Psi \qquad (A.5) $$

$$ \frac{z_{\star}}{R_{\text{orb}}} = - \sin(\Omega_{\text{orb}}t + \varphi_0) \sin \Psi. \qquad (A.6) $$

Figure 3 shows the trajectories of the planet GJ 436b as seen in the reference frame of the star during one rotational period of the star (about 44 days)

We note that, in the case of aligned systems, $\Psi = 0$ and we have

$$ x_{\star} = R_{\text{orb}} \cos(\Omega_{\text{orb}}t + \varphi_0 + \Omega_{\star}t) \qquad (A.7) $$

$$ y_{\star} = R_{\text{orb}} \sin(\Omega_{\text{orb}}t + \varphi_0 + \Omega_{\star}t) \qquad (A.8) $$

$$ z_{\star} = 0. \qquad (A.9) $$

Said differently, in aligned systems, the planet re-encounters the same stellar wind property once every $2\pi/(\Omega_{orb} + \Omega_{\star})$, as long as the stellar magnetic field does not evolve significantly.

## Appendix B: Interconnecting loop model

A third star-planet scenario we also explore in this paper is based on the idea presented by Lanza (2013), in which a magnetic loop connecting the star and the planet would suffer a rupture after being stretched by the planet. The idea is that one footprint of the loop sits at the surface of the star, while the other footprint lies at the surface of the planet. As the planet moves through its orbit, the magnetic loop is stretched (stressed) and when it breaks, it can release energy, which travels towards the star, giving rise for example to anomalous hot spots at the stellar surface (Shkolnik et al. 2005; Cauley et al. 2019). This is the scenario that Strugarek et al. (2022) named 'stretch and break'; in their most recent work on star-planet interactions in the GJ 436 system, this model was also investigated by Loyd et al. (2023). In this scenario, the model assumes that the stellar closed corona extends up to the distance of the planet, which is not the case in the stellar wind models we presented here – as shown in Sect. 3, the stellar wind stretches open the magnetic field lines of the star at distances within the planetary orbital distance. Additionally, the interconnecting loop scenario assumes that the stellar
## Appendix A: Calculation of orbital trajectories for a misaligned orbit

To calculate the orbital trajectories in the reference frame of the star, we compute two coordinate transformations: one from the orbital reference frame to the observer's frame, and another from the observer's frame to the stellar frame. These frames are illustrated in Fig. A.1.

We start from the orbital reference frame {$x_{pl}, y_{pl}, z_{pl}$}, such that the $x_{pl}y_{pl}$-plane contains the orbit of the planet and the orbital axis is along $z_{pl}$. The origin of this reference frame is centred on the star. In this reference frame, the (circular) planetary motion is simply described as

$$ \begin{bmatrix} x_{pl} \\ y_{pl} \\ z_{pl} \end{bmatrix} = \begin{bmatrix} R_{orb} \cos \varphi_{pl} \\ R_{orb} \sin \varphi_{pl} \\ 0 \end{bmatrix}, \qquad (A.1) $$

where the phase of the orbit is $\varphi_{pl} = \Omega_{orb}t + \varphi_0$, with $t$ describing the time and $\Omega_{orb} = 2\pi/P_{orb}$ the orbital rotation rate. Here, we assume a non-null initial phase $\varphi_0$, which we computed later in this section to find the location of GJ 436b at zero phase of the surface magnetic map, which corresponds to the first spectropolarimetric observation as reported in Bellotti et al. (2023).

Now we take the observer's (inertial) reference frame {$x_{obs}, y_{obs}, z_{obs}$}, whose origin is also centred on the star, and $z_{obs}$ is aligned with the stellar rotation. The observer is located at positive $x_{obs}$ and we assume that $x_{obs} \parallel x_{pl}$. The angle between the orbital spin and the stellar spin axis is the spin-orbit angle, $\Psi$; thus, the angle between $z_{obs}$ and $z_{pl}$ is also $\Psi$. In this reference frame, the orbital motion of the planet is described as

$$ \begin{bmatrix} x_{obs} \\ y_{obs} \\ z_{obs} \end{bmatrix} = \begin{bmatrix} 1 & 0 & 0 \\ 0 & \cos \Psi & \sin \Psi \\ 0 & -\sin \Psi & \cos \Psi \end{bmatrix} \begin{bmatrix} x_{pl} \\ y_{pl} \\ z_{pl} \end{bmatrix}, \qquad (A.2) $$

where we performed an anti-clockwise rotation around the axis $x_{obs} \equiv x_{pl}$ by an angle $\Psi$. According to our reference frame, the mid-transit occurs when the planet is at {$x_{obs}, y_{obs}, z_{obs}$} = {$R_{orb}, 0, 0$}.

<figure><img src="image_7.png" /><figcaption>Fig. A.1. Reference frames used in the coordinate transformation, with the star given by the orange circle, the planet indicated by the red circle, and its orbital path shown in red. The planet's reference frame {$x_{pl}, y_{pl}, z_{pl}$} (in red) contains the orbital plane of the planet $x_{pl}y_{pl}$ and orbital spin axis along $z_{pl}$. The observer's inertial reference frame {$x_{obs}, y_{obs}, z_{obs}$} (in blue) is centred on the star, with the observer located at $+x_{obs}$ and the stellar rotation axis aligned along $z_{obs}$. The reference frame {$x_{\star}, y_{\star}, z_{\star}$} (in black) co-rotates with the star. The spin-orbit angle is $\Psi$.</figcaption></figure>

The reference frame co-rotating with the star {$x_{\star}, y_{\star}, z_{\star}$} is also centred on the star, with the stellar spin axis along $z_{\star}$, which is parallel to $z_{obs}$. To transform from the observer's reference frame to the stellar co-rotating frame, we performed a clockwise rotation around the axis $z_{\star} \equiv z_{obs}$ by an angle $\Omega_{\star}t$:

$$ \begin{bmatrix} x_{\star} \\ y_{\star} \\ z_{\star} \end{bmatrix} = \begin{bmatrix} \cos(\Omega_{\star}t) & -\sin(\Omega_{\star}t) & 0 \\ \sin(\Omega_{\star}t) & \cos(\Omega_{\star}t) & 0 \\ 0 & 0 & 1 \end{bmatrix} \begin{bmatrix} x_{obs} \\ y_{obs} \\ z_{obs} \end{bmatrix}. \qquad (A.3) $$

Therefore, to derive the orbital path of the planet in the reference frame of the star, we substitute (A.1) and (A.2) into (A.3), thus obtaining

$$ \frac{x_{\star}}{R_{\text{orb}}} = \cos(\Omega_{\text{orb}}t + \varphi_0) \cos(\Omega_{\star}t) - \sin(\Omega_{\text{orb}}t + \varphi_0) \sin(\Omega_{\star}t) \cos \Psi \qquad (A.4) $$

$$ \frac{y_{\star}}{R_{\text{orb}}} = \cos(\Omega_{\text{orb}}t + \varphi_0) \sin(\Omega_{\star}t) + \sin(\Omega_{\text{orb}}t + \varphi_0) \cos(\Omega_{\star}t) \cos \Psi \qquad (A.5) $$

$$ \frac{z_{\star}}{R_{\text{orb}}} = -\sin(\Omega_{\text{orb}}t + \varphi_0) \sin \Psi. \qquad (A.6) $$

Figure 3 shows the trajectories of the planet GJ 436b as seen in the reference frame of the star during one rotational period of the star (about 44 days)

We note that, in the case of aligned systems, $\Psi = 0$ and we have

$$ x_{\star} = R_{\text{orb}} \cos(\Omega_{\text{orb}}t + \varphi_0 + \Omega_{\star}t) \qquad (A.7) $$

$$ y_{\star} = R_{\text{orb}} \sin(\Omega_{\text{orb}}t + \varphi_0 + \Omega_{\star}t) \qquad (A.8) $$

$$ z_{\star} = 0. \qquad (A.9) $$

Said differently, in aligned systems, the planet re-encounters the same stellar wind property once every $2\pi/(\Omega_{orb} + \Omega_{\star})$, as long as the stellar magnetic field does not evolve significantly.

## Appendix B: Interconnecting loop model

A third star-planet scenario we also explore in this paper is based on the idea presented by Lanza (2013), in which a magnetic loop connecting the star and the planet would suffer a rupture after being stretched by the planet. The idea is that one footprint of the loop sits at the surface of the star, while the other footprint lies at the surface of the planet. As the planet moves through its orbit, the magnetic loop is stretched (stressed) and when it breaks, it can release energy, which travels towards the star, giving rise for example to anomalous hot spots at the stellar surface (Shkolnik et al. 2005; Cauley et al. 2019). This is the scenario that Strugarek et al. (2022) named 'stretch and break'; in their most recent work on star-planet interactions in the GJ 436 system, this model was also investigated by Loyd et al. (2023). In this scenario, the model assumes that the stellar closed corona extends up to the distance of the planet, which is not the case in the stellar wind models we presented here – as shown in Sect. 3, the stellar wind stretches open the magnetic field lines of the star at distances within the planetary orbital distance. Additionally, the interconnecting loop scenario assumes that the stellar
<figure><img src="image_7.png" /><figcaption>**Fig. B.1.** Interconnecting loop model. Top: Magnetic field of the star at the orbital distance of the planet, assuming a potential field extrapolation, with a source surface of ≈ 15R<sub>*</sub> (in this case, no wind model is used). Bottom: Maximum dissipated power in the interconnecting loop model, assuming a polar planetary magnetic field strength (dipole) of 2 G (dashed) and 6 G (solid). The colour scheme follows the colour bar shown in Fig. 3. The grey arrows in the upper horizontal axes indicate approximate times of mid-transits.</figcaption></figure>

With this, we can then compute the power released in the interconnecting loop model using (Lanza 2013)

$$ \mathcal{P}_{\text{sb}} = (2\pi R_p^2) f_{\text{open}} \left[ \frac{B_p^2 v_K}{4\pi} \right], \qquad (\text{B.1}) $$

where the Poynting flux (term within brackets) is now across the planet. The fraction of the planetary surface area that has open magnetic field lines (Adams 2011) and thus have interconnecting magnetic loops is

$$ f_{\text{open}} = 1 - \left(1 - \frac{3\zeta^{1/3}}{2 + \zeta}\right)^{1/2}, \qquad (\text{B.2}) $$

with $\zeta = B/B_p$. Here, B is the stellar magnetic field extrapolated out to the orbit of GJ 436b and $B_p$, as before, is the assumed polar magnetic field strength of a dipolar planetary magnetic field. In the interconnecting loop scenario, the atmosphere of the planet is assumed to be ionised down to the surface, where the planetary magnetic field is $B_p$. If the atmosphere becomes neutral at a certain height above the surface, then the value of $B_p$ is that at such a height, which is smaller than the surface field. Therefore, Eq. (B.1) provides an upper limit of the power released in the interaction. Assuming $B_p = 2$ G, we find that the power released in the interconnecting loop model is on average $1.7 \times 10^{25}$ erg s<sup>-1</sup>. For a magnetic field of $B_p = 6$ G, the power released increases to $10 \times 10^{25}$ erg s<sup>-1</sup> (this value is similar to the maximum power derived by Loyd et al. 2023). This is shown in the bottom panel of Fig. B.1.

We confirm the findings of previous studies (Cauley et al. 2019; Strugarek et al. 2022), which showed that the interconnecting loop scenario provides the largest maximum power among the three scenarios explored in this paper (see the first two scenarios in Sect. 4), and is also in line with the estimates provided in Lanza (2013). According to the interconnecting loop scenario, a planetary magnetic field $B_p \le 6$ G can reproduce the maximum power observed in star–planet interactions in the GJ 436 system (Loyd et al. 2023).

magnetic field is close to a potential field, so that it can topologically reconnect with the magnetic field of the planet that is also potential. In the presence of a stellar wind though, the coronal magnetic field departs from a potential state (Vidotto et al. 2011b).

To explore the interconnecting loop scenario, we ran a model without the presence of a stellar wind, in which we assume that the stellar magnetic field remains potential and closed up to the orbit of the planet. To compute the field strength at the orbit of the planet, we use a potential field source surface extrapolation. In this model, the source surface represents the distance above which the magnetic field lines become open. To make sure the stellar magnetic field lines are closed up to the planet’s orbit, we chose a source surface of ≈ 15R<sub>*</sub> (i.e. beyond the planet’s orbit). Using the surface magnetic map from Fig. 1 and our potential field model, the magnetic field of the star at the orbit of the planet is shown in the top panel of Fig. B.1, where we see that the coronal magnetic field is on average $7.5 \times 10^{-3}$ G.
<figure><img src="image_8.png" /><figcaption>**Fig. B.1.** Interconnecting loop model. Top: Magnetic field of the star at the orbital distance of the planet, assuming a potential field extrapolation, with a source surface of ≈ 15R<sub>⋆</sub> (in this case, no wind model is used). Bottom: Maximum dissipated power in the interconnecting loop model, assuming a polar planetary magnetic field strength (dipole) of 2 G (dashed) and 6 G (solid). The colour scheme follows the colour bar shown in Fig. 3. The grey arrows in the upper horizontal axes indicate approximate times of mid-transits.</figcaption></figure>

With this, we can then compute the power released in the interconnecting loop model using (Lanza 2013)

$$ \mathcal{P}_{\text{sb}} = (2\pi R_p^2) f_{\text{open}} \left[ \frac{B_p^2 v_K}{4\pi} \right], \qquad (\text{B.1}) $$

where the Poynting flux (term within brackets) is now across the planet. The fraction of the planetary surface area that has open magnetic field lines (Adams 2011) and thus have interconnecting magnetic loops is

$$ f_{\text{open}} = 1 - \left(1 - \frac{3\zeta^{1/3}}{2 + \zeta}\right)^{1/2}, \qquad (\text{B.2}) $$

with $\zeta = B/B_p$. Here, B is the stellar magnetic field extrapolated out to the orbit of GJ 436b and $B_p$, as before, is the assumed polar magnetic field strength of a dipolar planetary magnetic field. In the interconnecting loop scenario, the atmosphere of the planet is assumed to be ionised down to the surface, where the planetary magnetic field is $B_p$. If the atmosphere becomes neutral at a certain height above the surface, then the value of $B_p$ is that at such a height, which is smaller than the surface field. Therefore, Eq. (B.1) provides an upper limit of the power released in the interaction. Assuming $B_p = 2$ G, we find that the power released in the interconnecting loop model is on average $1.7 \times 10^{25}$ erg s<sup>-1</sup>. For a magnetic field of $B_p = 6$ G, the power released increases to $10 \times 10^{25}$ erg s<sup>-1</sup> (this value is similar to the maximum power derived by Loyd et al. 2023). This is shown in the bottom panel of Fig. B.1.

We confirm the findings of previous studies (Cauley et al. 2019; Strugarek et al. 2022), which showed that the interconnecting loop scenario provides the largest maximum power among the three scenarios explored in this paper (see the first two scenarios in Sect. 4), and is also in line with the estimates provided in Lanza (2013). According to the interconnecting loop scenario, a planetary magnetic field $B_p \le 6$ G can reproduce the maximum power observed in star–planet interactions in the GJ 436 system (Loyd et al. 2023).

magnetic field is close to a potential field, so that it can topologically reconnect with the magnetic field of the planet that is also potential. In the presence of a stellar wind though, the coronal magnetic field departs from a potential state (Vidotto et al. 2011b).

To explore the interconnecting loop scenario, we ran a model without the presence of a stellar wind, in which we assume that the stellar magnetic field remains potential and closed up to the orbit of the planet. To compute the field strength at the orbit of the planet, we use a potential field source surface extrapolation. In this model, the source surface represents the distance above which the magnetic field lines become open. To make sure the stellar magnetic field lines are closed up to the planet's orbit, we chose a source surface of ≈ 15R<sub>⋆</sub> (i.e. beyond the planet's orbit). Using the surface magnetic map from Fig. 1 and our potential field model, the magnetic field of the star at the orbit of the planet is shown in the top panel of Fig. B.1, where we see that the coronal magnetic field is on average $7.5 \times 10^{-3}$ G.