Topology and obliquity of core magnetic fields in shaping seismic properties of slowly rotating evolved stars

It is thought that magnetic fields must be present in the interiors of stars to resolve certain discrepancies between theory and observation (e.g. angular momentum transport), but such fields are difficult to detect and characterise. Asteroseismology is a powerful technique for inferring the internal structures of stars by measuring their oscillation frequencies, and succeeds particularly with evolved stars, owing to their mixed modes, which are sensitive to the deep interior. The goal of this work is to present a phenomenological study of the combined effects of rotation and magnetism in evolved stars, where both are assumed weak enough that first-order perturbation theory applies, and we focus on the regime where Coriolis and Lorentz forces are comparable. Axisymmetric"twisted-torus"field configurations are used, which are confined to the core and allowed to be misaligned with respect to the rotation axis. Factors such as the field radius, topology and obliquity are examined. We observe that fields with finer-scale radial structure and/or smaller radial extent produce smaller contributions to the frequency shift. The interplay of rotation and magnetism is shown to be complex: we demonstrate that it is possible for nearly symmetric multiplets of apparently low multiplicity to arise even under a substantial field, which might falsely appear to rule out its presence. Our results suggest that proper modelling of rotation and magnetism, in a simultaneous fashion, may be required to draw robust conclusions about the existence/non-existence of a core magnetic field in any given object.


INTRODUCTION
Magnetism in stars can occur in regions where a dynamo is operating to actively generate the field, or they may be fossil fields, i.e. accreted passively from the parent gas cloud (Mestel 2012). Regions in a star where dynamo action can operate are those unstable to turbulent convection, which occurs in the cores of intermediateto high-mass (> 1.2 ⊙ ) main sequence stars, the envelopes of less massive main sequence stars, and the envelopes of red giants of all masses (Maeder 2008). While it is possible to directly observe magnetic fields that penetrate the surface (e.g. through spectropolarimetry), whether this be from dynamo action or a fossil field (Donati & Landstreet 2009;Wade et al. 2016;Shultz et al. 2018), core magnetism is much harder to establish and probe. For this, indirect means have been sought in the form of numerical simulations (Brun et al. 2005;Featherstone et al. 2009), and asteroseismology (Fuller et al. 2015;Stello et al. 2016).
Asteroseismology is the technique of deducing a star's interior properties from analysis of its natural modes of oscillation. Applied first to the Sun (helioseismology), it has since been turned with great success to more distant stars, yielding numerous breakthroughs par-★ E-mail: stl36@cam.ac.uk ticularly in the context of evolved stars (subgiants and red giants). This owes to the existence of mixed modes in such objects, which involve a coupling of fluid motions deep in the core (gravity waves/gmode oscillations) to those at the surface (acoustic waves/p-mode oscillations). This has enabled core rotation rates (Beck et al. 2012) and the presence of core helium burning (Bedding et al. 2011) to be inferred. Signatures of strong (i.e. dynamically significant) magnetic fields in the cores of evolved stars are also thought to manifest in the form of mode depression (Mosser et al. 2012a), where resonant interactions between gravity and Alfvén waves cause mode conversion and energy dissipation by phase mixing (Lecoanet et al. 2017;Loi & Papaloizou 2018;Loi 2020a). However, the exact origin of observed mode depression and its necessity to invoke magnetic fields remains controversial (García et al. 2014;Mosser et al. 2017).
Numerous open questions surround the existence and role of internal magnetic fields in stellar physics/evolution. It is speculated that they may be able to explain observed efficiencies of angular momentum transport ), which are much higher than non-magnetic processes can collectively account for, although the exact mechanism(s) are under debate (Fuller et al. 2019;den Hartogh et al. 2020;Takahashi & Langer 2020). Magnetic fields may also help resolve the problem of photospheric chem-ical abundances in evolved stars, where enhanced mixing within radiative zones is required (Mathis & Zahn 2005;Busso et al. 2007). It is therefore of interest to find ways of detecting and characterising magnetic fields in deep stellar interiors, for which mixed-mode asteroseismology of evolved stars offers great promise.
The mode depression phenomenon, which is exclusive to those red giant stars >1.2 ⊙ (previously able to host core dynamos when on the main sequence), occurs in only a fraction of these (∼50% for those >1.6 ⊙ ; Stello et al. 2016). It is therefore to be wondered what explains the remaining stars: might they also have significant core fields but just below the critical threshold needed for mode depression, or is there a genuine dichotomy in the field strengths/properties? Since magnetic fields are not spherically symmetric, the Lorentz force lifts the degeneracy of modes with the same radial order and spherical harmonic degree ℓ but different azimuthal order , producing frequency splitting and giving rise to a multiplet. Similarly, the Coriolis force (arising from rotation) also produces frequency splitting (Beck et al. 2012;Mosser et al. 2015). At low rotation rates and/or weak fields, this problem can be treated using first-order perturbation theory, which calculates the frequency correction to the associated unperturbed eigenmode assuming that the Coriolis and Lorentz forces are small compared to pressure and buoyancy. This approach has been used to treat solar p-modes (Gough & Thompson 1990), fundamental/low-order p-modes in Cepheid variables (Shibahashi & Aerts 2000), and rmodes in degenerate stars (Morsink & Rezania 2002) under the influence of rotation and magnetism. For g-modes, the effects of a magnetic field were investigated by Hasan et al. (2005) for slowly pulsating B-stars, and Rashba et al. (2007) in the case of the Sun.
At larger rotation rates/field strengths, e.g. where the rotation/Alfvén frequency is a substantial fraction of the mode frequency, non-perturbative approaches must be used. For example, there is a body of recent work on rapidly rotating stars where the Coriolis force is included via the traditional approximation (Buysschaert et al. 2018;Prat et al. 2019Prat et al. , 2020, and first-order perturbation theory is used for the Lorentz force only. This is relevant for massive main sequence stars, for which magnetic braking does not efficiently operate to slow down their rotation. In contrast, evolved stars rotate much more slowly due to the huge expansion of their envelopes upon leaving the main sequence, and so rotational effects may be treated perturbatively. As mentioned above, the existence/properties of magnetic fields in red giant cores is still largely speculative. While a sizeable fraction may have strong fields that explain their mode depression, detailed predictions of how weaker magnetism might impact the seismic properties of evolved stars have been lacking. In particular, there are few studies examining the simultaneous consequences of rotation and magnetism in these objects, although several works have attempted to quantify the effects of magnetism by itself Gomes & Lopes 2020;Loi 2020b). The recent works of Mathis et al. (2021) and Bugnet et al. (2021), which have included the effects of aligned rotation in addition to a large-scale axisymmetric field, have begun to fill these gaps. However, all the above works have assumed field configurations that are simple and large scale. In contrast, the current work aims to treat the unexplored regime of more complex radial field structure, that may also be misaligned with respect to rotation (i.e. non-axisymmetric).
It is the goal of this work to perform a phenomenological study into the combined effects of rotation and magnetism on mixed modes in evolved stars, where both Coriolis and Lorentz forces are assumed to be weak enough that first-order perturbation theory applies. The focus lies on low-degree mixed modes of short radial wavelength, and the regime where the two forces are comparable in strength. Axisymmetric "twisted torus" configurations for the magnetic field will be used, and allowed to be inclined with respect to the rotation axis. This is motivated by the knowledge that rotation and magnetic axes often do not coincide, as seen both in stars with fossil fields (Henrichs et al. 2013;Braithwaite & Spruit 2017) and in objects with active dynamos including the Sun (Gosling 2007) and the Earth (Merrill 2010). Rotation will be treated as uniform for simplicity; note that in the core region this would be physically justified as magnetic fields tend to enforce solid-body rotation (Spruit 1999). Besides obliquity of the field, the impact of topology will also be examined, by considering twisted torus configurations of more complex radial structure. Note that no good knowledge exists about the most likely field topologies that may be found in red giant cores. Speculatively, given that the most well-observed dynamos exhibit periodic reversals on relatively short timescales (Jacobs 1994;Ossendrĳver 2003), and simulations of stellar core dynamos show similar behaviour (Brun et al. 2005), one might imagine that as the convective core recedes over the main sequence it might leave behind magnetised shells of alternating sign. Whether or how such small-scale complex structure (if it can survive to later stages) might influence frequency splittings has not previously been studied. It is to be noted that configurations with larger spatial scales are lower in energy and expected to be more stable (Broderick & Narayan 2008;Duez & Mathis 2010). However, we consider here the possibility that some dynamos may preferentially create fields of smaller-scale structure, which may or may not eventually collapse to a lower energy state. This paper is structured as follows. In Section 2, we introduce the stellar models and magnetic field configurations. In Section 3, we describe how the basic (unperturbed) eigenmodes were computed, and review relevant aspects of first-order perturbation theory. Results are presented in Section 4 and discussed further in Section 5, which covers implications for asteroseismic inference and limitations of the framework. Finally, we conclude in Section 6.

MODELS
Four stellar models were examined in this study, all of mass * = 2 ⊙ . Two are polytropes of differing index, and two are realistic evolved stellar models generated by the publicly available stellar evolutionary code 'Modules for Experiments in Stellar Astrophysics' ( , version r11701) (Paxton et al. 2011). Their parameters are summarised in Table 1, and described in further detail below. The structure of these models was computed neglecting rotation and magnetism, as it is assumed that these effects are too small to produce significant deformations away from sphericity, or otherwise influence the hydrostatic background structure. The magnetic field models were calculated separately (see Section 2.3 for details), scaled to a desired strength and imposed on the core region of each model within a boundary of radius f .
The rotation profile is assumed to be uniform and described by a single scalar frequency Ω, whose value is indicated in Table  1 for each of the models. For the polytropic models (A and B), Ω was chosen to be a small fixed fraction (0.2%) of the dynamical frequency dyn = * / 3 * , where * is the stellar radius. For the models (C and D), Ω was chosen to be consistent with characteristic rotational splittings of ∼100 nHz observed for dipole modes in red giant stars of around 2 ⊙ (Mosser et al. 2012b, fig. 6). This translates to Ω values of about 1% and 3% of the dynamical frequencies of Models C and D, respectively.

Polytropes
Solutions to the Lane-Emden equation which comes from substituting the polytropic relation ∝ 1+1/ into the equations of hydrostatic equilibrium, provide simple models for stars. Here is the gas pressure, is the mass density, is the polytropic index, is the radial coordinate, and is the polytropic temperature. These models are easy to generate and were considered in this study to see how much the results depended on the structure of the background. Two polytropic models were generated, having the same mass * = 2 ⊙ and radius * = 6 ⊙ , but different values of . These are referred to as Models A and B, where A has = 4.2 and B has = 4.6. The value of controls the central condensation of a polytrope: Model B has a smaller, denser core and larger, more tenuous envelope than Model A. These are shown in Fig. 1, which shows the radial variation of , and the Lamb and buoyancy frequencies, ℓ and . A key difference between polytropes and realistic stellar models is that the former lack an evanescent zone between regions of p-mode and g-mode propagation, and are nonconvective throughout. Thus all modes in polytropes are heavily mixed, with substantial amounts of both p-and g-like character.

Evolved stellar models
A sequence of 2 ⊙ stellar models (with metallicity = 0.02) along a single evolutionary track was generated using , the inlist of which is given in Appendix A. Two of the output models were selected for further analysis, corresponding to a subgiant (Model C) and a red giant (Model D), at ages of 976 Myr and 1.01 Gyr, respectively. It was noticed that under the default settings of a significant level of jaggedness on the grid scale would be present. To remedy this, the mesh criteria were tightened to force to use a finer grid, and further smoothing was applied to the outputs. The resultant profiles of , and for Models C and D are shown in Fig. 2.
In more detail, the smoothing procedure involved first replacing each point by that lying on the least-squares line of best fit through the 13 neighbouring points (cf. boxcar smoothing), which were then linearly interpolated to a uniform grid having the same number of points. This was downsampled by a factor of 5, and then spline interpolated to a uniform grid having a chosen final number of points, which were 10 5 and 10 6 points for Models C and D, respectively. Note that the values listed under '# grid points' in Table  1 are only half as large, as they refer to the size of the grid on which the eigenmodes were saved (this stems from the way the integration routine was implemented).

Magnetic field
Realistic modelling of magnetic fields in stars is a non-trivial task, owing to the need to satisfy various physical and stability constraints. Besides the well-known solenoidal condition ∇ · B = 0, where B is the magnetic field, all components of B need to be finite and continuous everywhere (to avoid infinite current sheets). If one desires to model a spatially confined field (e.g. core of an evolved star) then such a field cannot in general be force-free (Spruit 2013). Furthermore, any purely poloidal or purely toroidal configuration is unstable and so only mixed configurations are allowed (Tayler 1973;Flowers & Ruderman 1977). Most simple configurations, including uniform vertical fields, purely toroidal fields, and dipole fields, which have been widely used in many works, violate one or more of the above and are thus unlikely to be good descriptions of fields in reality. Instead, various analytical studies and numerical simulations point strongly towards so-called "twisted-torus" configurations (Prendergast 1956;Braithwaite & Nordlund 2006;Yoshida et al. 2006;Duez & Mathis 2010;Duez et al. 2010b), which are axisymmetric and dipole-like in angular appearance (in the sense that they have two poles) but differ in several important ways from actual dipoles, namely that they have no central singularity, vanish smoothly at a finite radius, and possess a stabilising toroidal component of comparable magnitude to the poloidal component.
Such configurations, despite being moderately complex in appearance, are surprisingly easy to construct analytically. This was first achieved by Prendergast (1956) for incompressible fluids and later generalised to the compressible case by Duez & Mathis (2010) and Duez et al. (2010b). The derivation will not be covered here; the reader is encouraged to refer to the above works for details, also see Loi (2020a) for a summary. Although derived in the non-rotating case, numerical work by Duez (2011) finds similar configurations in the presence of rotation, and these to be preferentially oblique with respect to the rotation axis. Such "twisted-torus" field configurations will be referred to here as Prendergast solutions.
Each Prendergast solution is completely described by two quantities: a scalar function of radius, Ψ( ), known as the radial flux function, and a parameter . These yield the magnetic field components according to where ( , ′ , ′ ) are spherical polar coordinates defined with respect to the magnetic axis. Primes are used to distinguish the angular coordinates from those defined with respect to the rotation axis (i.e. and ), which is allowed to be different from the magnetic axis. The parameter must be a root of where 1 is a spherical Bessel function of the first kind. Roots of (4) will exist if does not decrease too rapidly compared to f , and in general there can be multiple roots due to the oscillatory nature of 1 . Roughly speaking, gives the inverse length scale of Ψ, and so larger roots correspond to configurations with progressively finer radial structure. Once is chosen, Ψ is constructed as where ( 1 , 2 ) ≡ 1 ( 2 ) 1 ( 1 ) − 1 ( 1 ) 1 ( 2 ), and 1 is a spherical Bessel function of the second kind. The final desired field strength determines the scaling of Ψ.
Prendergast solutions for the first nine roots of Model C are shown in Fig. 3. While the lowest-order configuration may be described as a twisted torus, higher-order configurations take the form of multiple nested tori. Note how the north-south projection of the field lines flips sign between radial shells, which one might imagine could mimic the end result of a periodically reversing dynamo. Similar plots for the other stellar models, and corresponding distributions of the total field strength, are included as Supplementary  Fig. 1, but for the two models (C and D). Note that plots of all quantities are zoomed in to the core region; the full stellar radii are 4.3 and 7.8 ⊙ for Models C and D, respectively.
Figs S1-S7. For the most part, the field strength is maximal at the centre, with additional smaller local maxima at larger radii. Along the boundary = f , all field components smoothly match to the zero solution.
The default values of f for each of the models are listed in Table 1, and marked in Figs 1 and 2. In the case of the models, they were chosen to lie within the boundary of the former convective core, which was identified through inspection of the H and He composition profiles. The overall scaling of the field strength was controlled through the central Alfvén speed, the default values of which are also listed in Table 1. These were chosen such that the magnetic contribution to the frequency splitting would be about the Table 1. Summary of parameters for the four stellar models, including their mass * , radius * , dynamical frequency dyn = * / 3 * , dynamical speed dyn = * / * , grid size, range of oscillation frequencies considered, number of modes found, rotation frequency Ω, radial extent of the magnetic field f , central Alfvén speed A,cen , corresponding value of the central field strength cen , and an estimate of the critical strength crit for comparison. Note that values of f and A,cen shown are defaults. In some parts of this paper differing values are used, but where this occurs it will be explicitly stated. Otherwise, the default values should be assumed.

Basic eigenmodes
First, we need to obtain the set of basic eigenmodes that exist in the absence of rotation and magnetism. The time-dependent fluid displacement vector field for a normal mode of frequency can be written where Note that torsional motions have been neglected as the only restoring forces at this stage are pressure and buoyancy. Also, for spherically symmetric backgrounds ℓ, and ℓ, will not depend on , so we shall write ℓ, → ℓ , ℓ, → ℓ . Substituting into the equations of motion and linearising, and projecting out a single spherical harmonic, we get the equations of Note that we have made the Cowling approximation (suitable for modes with small radial scales) thus reducing the system to second order, and neglected spatial variations of the adiabatic index . The eigenmodes used for subsequent analysis in this work were obtained by directly solving Equations 8 and 9 under appropriate boundary conditions (regularity at = 0, and vanishing Lagrangian pressure perturbation at the surface). They were solved via standard methods (shooting with a fourth-order Runge-Kutta scheme and matching at an intermediate radius) to obtain a set of unperturbed eigenmodes for spherical harmonic degrees ℓ = 1, 2 and 3 over chosen ranges. These were [0.5, 1] dyn for Models A and B, and [8, 12] dyn for Models C and D. The numbers of modes found within each interval are listed in Table 1. The radial order of the various modes was computed using the Eckart (1960) scheme. Figure 4 plots the frequencies and radial orders of all modes. The vast majority of these are highly g-dominated, with large negative radial orders and frequencies that scale roughly as 1/| |. Examples of the associated eigenfunctions are shown in Fig. 5, where the finest scales of variation can be seen to occur in the deep interior. A useful measure of the p-or g-like character of a mode is the mode inertia which measures the mass of fluid displaced and is larger for modes localised to the core where is larger, i.e. g-dominated modes. Figure 6 shows the values of ℓ for all modes of Models C and D. The small number of modes with very low inertia correspond to the p-dominated modes, which occur near the frequencies of pure envelope p-modes. For Models A and B, their lack of an evanescent zone means that ℓ values are relatively constant over the considered frequency range; these are not shown.

First-order perturbation theory
To calculate frequency corrections to the basic eigenmodes induced by rotation and magnetism, first-order perturbation theory was used, which is a standard formalism e.g. see Dziembowski & Goode (1984); Unno et al. (1989). It is valid in the limit of small perturbing forces, which is the regime focused on here. This section overviews relevant aspects of this theory, in particular its application to the problem of misaligned rotation and magnetism. The presence of a non-spherically symmetric force (e.g. Coriolis, Lorentz) lifts the degeneracy otherwise possessed by modes of the same and ℓ but different . For a given ℓ, there exist 2ℓ + 1 values of going from −ℓ to +ℓ. To first order, pure rotation of frequency Ω perturbs the mode frequencies by a value proportional to Ω, giving a multiplet of 2ℓ + 1 peaks. On the other hand, pure magnetism of characteristic field strength perturbs the frequencies by a value proportional to 2 and independent of the sign of , giving a multiplet of ℓ + 1 peaks where modes of the same | | continue to experience a degeneracy. Note that this statement applies to axisymmetric fields, for which the quantum numbers ℓ and can be defined, and not for general non-axisymmetric fields where ℓ and would lose their meaning.  Under the combined influence of the two effects, the picture is more complicated. From the point of view of the star (i.e. in the stationary/corotating frame), there can still be only 2ℓ+1 frequencies since that is the number of originally degenerate modes. Under axisymmetric conditions (aligned rotation and magnetic axes) each of these is associated with a different , and from the point of view of an observer (i.e. in the inertial frame) there are thus 2ℓ + 1 peaks, where the frequencies are additionally Doppler shifted by a value of Ω. This is similar to the pure rotation case, except that due to the Lorentz force the multiplet may exhibit asymmetries; for a complete discussion and seismic diagnosis of these asymmetries, we refer the reader to the detailed discussion in Bugnet et al. (2021). Under nonaxisymmetric conditions (misaligned rotation and magnetic axes), there would still be 2ℓ + 1 modes in the stationary/corotating frame, but now each mode can no longer be associated with a single since the system has no axis of symmetry. Rather, each would be an admixture of the 2ℓ + 1 different values of , giving rise to (2ℓ + 1) 2 peaks in the inertial frame after Doppler shifting. Their amplitudes are directly related to the coefficients of expansion: for example, in the limit of mutual alignment (which is a special case of this more general framework), a single dominates strongly for each corotating-frame mode and thus only 2ℓ + 1 out of (2ℓ + 1) 2 possible peaks have non-zero amplitudes. The mathematical treatment is as follows. Let and ′ be defined as the azimuthal quantum numbers with respect to the rotation and magnetic axes, respectively, and let be the obliquity angle. Let and¯ be the frequencies in the inertial and corotating frames, respectively, and thus related by¯ = + Ω (this definition means that < 0 modes are prograde). The equation of motion can be written where L rot and L mag , which correspond to the Coriolis and Lorentz forces, are assumed to be much smaller than L 0 , which represents the combination of pressure and buoyancy. Their functional forms in terms of the fluid displacement and background quantities are The goal is to solve for the first-order frequency corrections¯ 1 , where¯ = 0 +¯ 1 + · · · . Here 0 is the frequency associated with the unperturbed eigenmode 0 , where = 0 + 1 +· · · . Substituting these expansions into (11) and truncating at first order gives where L 1 = L rot + L mag . Now L 0 is self-adjoint, implying that { 0 } form a complete orthogonal basis. Therefore, it is possible to express any 1 = ( ) 0 , where is some enumeration of the set of { 0 }. Combining this with the fact that 2 0 0 = L 0 0 , for some eigenmode labelled (with chosen and ℓ).
If L 1 were axisymmetric, it would suffice to consider ( ) 0 as being made up of a single . Then, we would be able to take the inner product of (16) with ( ) 0, , which is defined to be (6) but restricted to a fixed , and use orthogonality to get The inner product is defined as , = ∫ * · d , where the integral is over the volume of the star. However, if ≠ 0 then L 1 will not be axisymmetric, even if L rot and L mag individually are. In this case L 1 mixes the different , and (17) needs to be replaced by a (2ℓ + 1) × (2ℓ + 1) matrix equation. Let us write  Fig. 4). Black and red correspond to the radial and horizontal components of the fluid displacement, ℓ and ℓ . The associated model, and ℓ are indicated in the panel headings. The mixed character of the modes is evident from the simultaneously substantial core and surface displacements in both cases. In the lower panel, which is for a more centrally condensed model, the inset shows a zoom-in to the core region so that the g-mode oscillations can be better seen.
where ′ 0, ′ has the same form as 0, but with ( , , ) → ( ′ , ′ , ′ ), and we have dropped the labels with the understanding that all quantities now pertain to the mode .
To transform between coordinate frames ( , , ) and ( , ′ , ′ ) we invoke the Wigner -matrix, denoted here by D. Its entries are given by where the summation is such that the factorials are non-negative. These satisfy which when substituted into (7) leads to  (10), for the two models. These show minima corresponding to where the mode becomes more p-dominated. In the more evolved model (bottom), the pure g-mode spectrum is much denser than the pure p-mode spectrum, leading to fewer minima in the mixed-mode spectrum per given number of modes, compared to Model C (top).
Hence (18) can be expressed as a = Da ′ , where a = ( : = −ℓ, · · · , +ℓ) and a ′ = ( ′ ′ : ′ = −ℓ, · · · , +ℓ). We can then recast (17) in the more general form of a matrix eigenvalue problem: where the Wigner -matrix has the property that D −1 = D ⊤ , and It is easy to see that this reduces to the aligned case when = 0 =⇒ D = I, in which case the eigenvectors a are simply the columns of I. The vectorial nature of a is not to be confused with a spatial coordinate vector; rather, it is a coordinate representation of the admixture of , with 2ℓ + 1 components.
A word of caution is to be made about the self-adjointness of the Lorentz operator L mag . In the fully self-consistent case where deformation is taken into account, the combined operator in (11) would be self-adjoint thus producing real frequencies, which is expected in ideal MHD. Note that in practice this deformation is small, and its dynamical consequences can largely be neglected (Duez et al. 2010a). As it stands, L 0 and L rot defined in (12) and (13) are self-adjoint, but L mag as defined in (14) is not. The approximations used thus have the side effect of introducing imaginary components to the frequencies which are not physical. However, there exists a convenient remedy in the form of retaining only the first term on the RHS of (14), which is the dominant term for modes of short wavelength and also happens to be self-adjoint. Hence we make the approximation Also, since rotation frequencies are assumed to be small, we neglect the second (centrifugal) term in L rot in favour of the first (Coriolis): The L 0 term is left as is. With these simplifications, the elements of M rot and M mag are then where , , in (27) refer to the components of 0, , and ′ , ′ ′ , ′ ′ in (28) refer to the components of ′ 0, ′ . Here J ′ = ( ′ , ′ ′ , ′ ′ ) is the Eulerian perturbation to the current density. Expressions for the components of J ′ for the special case of the Prendergast solution can be found in Appendix B. Figures 7 and 8 show the inertial-frame frequency shifts 1 =¯ 1 − Ω for Models C and D as a function of 0 , for three different cases: pure rotation (black), rotation plus magnetic field where = 0 (red), and rotation plus magnetic field where = /4 (blue). Comparing the red curve with the black curve, two aspects can be noted: (i) pdominated multiplets are more symmetric than g-dominated ones, and (ii) the asymmetries are more pronounced at lower frequencies.

Comparison with pure rotation
The first point can be explained by the fact that the magnetic field is confined to the core where g-dominated mixed modes are localised; p-dominated mixed modes are thus more similar to modes in the case of pure rotation, for which 1 is roughly constant with 0 . The second point can be understood through the frequency dependence of the elements of M rot and M mag . While M mag ∝ 1/ 0 , the factors of 0 in M rot cancel out in the numerator and denominator, leaving no explicit 0 dependence. The increased importance of magnetic effects at lower frequencies can be physically understood in terms of the increase of the Alfvén frequency ∝ where is the wavenumber, thus bringing it closer to the mode frequency. Note that we are in the regime where ≪ 0 , and that ∼ 0 corresponds to the strong-field (dynamically significant) regime in which perturbation theory breaks down.
It is curious to note that for ℓ = 1, in both Models C and D, the blue curves ( = /4) exhibit more symmetric splittings than the red ( = 0). This appears to roughly hold across the whole frequency range. However, the ℓ = 2 case does not show this behaviour, at least for the value of tested.

Effect of obliquity
As mentioned in Section 3.2, an important consequence of obliquity is to mix the different components, giving rise to (2ℓ + 1) 2  Table 1) with = 0, and blue crosses are for the same field strength but with = /4. In all three cases the default rotation profile has been applied, which is Ω = 0.01 dyn for this model. The 2ℓ + 1 curves of the same colour in each panel correspond to the different modes in the multiplet; for the blue case the component with the largest | | has been selected for plotting. frequencies in the inertial frame. However, only 2ℓ + 1 blue curves have been chosen for plotting in Figs 7 and 8; these correspond to the components with the largest | | values, i.e. largest observed amplitudes. As the field strength and/or obliquity increase, it becomes less easy to identify the component with the largest | |, as the various | | values become comparable. The -mixing process is illustrated in Figs 9-11, for a chosen mode of Model C with ℓ = 1, 2 and 3. The values of Ω, A,cen and f are the same in all cases. Along the top row of each plot are the matrices M mag and M rot , whose entries quantitatively represent the corotating-frame frequency shift that would be induced by each effect in the absence of the other. Along the bottom row are matrices containing the 2ℓ + 1 eigenvectors of (22), where each column corresponds to one eigenvector. These are shown for two values of . When = 0 (aligned case) these yield the identity matrix (no -mixing), but for = /3 (oblique case) the off-diagonal components are non-zero. Comparing the three values of ℓ, it can be seen that the magnetic contribution to the frequency shifts, given by the components of M mag , increases for higher ℓ. Physically this can be understood from the fact that for modes possessing glike character, larger ℓ are associated with smaller spatial scales and therefore larger Alfvén frequencies, for the same 0 . In contrast, the rotational contribution, given by the components of M rot , decreases for increasing ℓ given the same . Consequently the amount ofmixing increases with ℓ, for the same obliquity. Figures 12-14 illustrate the splitting process for a chosen mode of Model C, for six different values of going from 0 (aligned) to /2 (perpendicular). Similar plots for Models A and D can be found in Supplementary Figs S8-S13. It can be seen that greater -mixing occurs for intermediate obliquities ( ∼ /4), higher ℓ and larger field strengths. In Fig. 12, it is apparent that in the inertial frame only 2ℓ + 1 out of (2ℓ + 1) 2 peaks dominate, consistent with the a matrix remaining approximately diagonal even at substantial obliquities (see Fig. 9). Closer inspection of how the shape of this dominant sub-multiplet changes with obliquity reveals that while this may be asymmetric in general, it is still possible for this to appear nearly symmetric for some values of and have a spacing close to the value expected of pure rotation (Fig. 12, lower left panel). However, the centroid of this sub-multiplet is offset compared to the pure rotational multiplet. At larger field strengths (see Fig. 14), there may be many peaks of significant amplitude in the inertial frame. Notably, due to the heavy mixing, it would be possible to find symmetric sub-multiplets among these (e.g. bottom row, middle panel).
Frequency splittings for all ℓ = 1 modes of Models C and D are summarised in Fig. 15, where points are coloured according to the associated value of | |. The increased -mixing with larger obliquity and field strength can clearly be seen across all modes, both p-and g-dominated. The main difference between p-and gdominated mixed modes appears to lie in the centroid frequency of the multiplet: this is located towards systematically higher values compared to the unperturbed value for g-dominated multiplets, a feature not present in the case of pure rotation (as seen in Figs 7 and 8), where the centroid offset is zero for all multiplets. This offset between p-and g-dominated multiplets increases towards lower frequencies, where magnetic effects become more significant.

Radial extent of the field
When constructing each Prendergast field solution, one has free choice over the parameter f , the radial extent of the field. The default values of f for each model used in preceding sections are  Table 1.  Figure 10. As for Fig. 9, but a mode with = −28, ℓ = 2. indicated in Table 1, but in this particular section we choose to vary f between 0.6-1.4 times the default, while keeping all other parameters (including A,cen ) the same. Figure 16 shows that the components of M mag , representing the magnetic contribution to the overall frequency shift, increase in rough proportion with f . This may be explained by the integral in (14) being larger when there is more volume of field to integrate over. However, this trend is not strict, and in some cases as for Model D the behaviour can be somewhat unpredicable. This appears to be tied to the complicating influence of variations in geometry/topology of the field as f is modified; this occurs for the Prendergast model because of the dependence of (4) and (5) on , whose shape over the interval ∈ [0, f ] does not scale straightforwardly with f . At points where the configuration acquires additional radial structure (see inset panels), this tends to disrupt the trend. Further discussion of the field topology can be found in the next section.
This plot also sheds light on the effect of central condensation. Consider Models A and B, which are both polytropic models with the same * and * (and therefore dynamical frequency), but have different central condensations due to their differing polytropic indices. They also have similar central field strengths (differing by a factor of two, being larger in Model A). However, the magnetic shifts in Model B are larger than A by a factor of about two, and since this scales with 2 , if Model B had the same central field strength then the frequency shift would be ∼8 times larger than in Model A. The increased importance of magnetic effects at higher central condensations can be attributed to the smaller spatial scales of oscillation in the core due to the larger buoyancy frequency, leading to a higher Alfvén frequency. This effect is more difficult to disentangle using Models C and D because the various modes differ wildly in their level of p-like versus g-like character, which represents a conflating influence.
A minor further comment about Fig. 16 concerns the ℓdependence for the different models. In Models A, B and C, for which the modes chosen for plotting all have relatively large amounts of mixed p-and g-like character, higher ℓ undergo larger shifts since the dominant factor is the shrinkage of spatial scales as ℓ increases (driving up the Alfvén frequency). However, in Model D, which is more evolved than Model C, the much greater p-dominated character for higher ℓ means that the behaviour is reversed: higher ℓ values undergo smaller magnetic shifts.
The inertial-frame frequency shifts (including the Coriolis force) for all modes of Model C are shown in Fig. 17, for two values of f . The g-dominated modes are seen to be much more sensitive to increases in f compared to p-dominated modes, as can be seen by their larger values of 1 . This can be attributed to their preferential localisation to the core, where the field is located.

Role of field topology
As mentioned in Section 2.3, multiple values of can satisfy (4) for a given combination of f and ( ). These correspond to fields of different topologies (cf. Fig. 3). Figure 18 shows how the components of M mag vary as a function of . The general trend is for these to decrease with increasing , suggesting that fields with more radial structure may produce smaller magnetic shifts, even if the central field strength is held constant. For Model D the trend is not so smooth (second and fourth values are outliers); inspection of Supplementary Figs S3 and S7 shows that these have unusual field structure, possessing a node/minimum very near the centre. The approach used here to control the field amplitude via scaling the central value thus produces anomalously large maximum field strengths and therefore frequency shifts. However, in all other cases the central and maximum field strengths are virtually identical.
Further insight into the cause of this systematic behaviour with comes from inspecting the spatial contribution of the integrand defined in (14), by considering the value of the integral evaluated from 0 to some finite rather than all the way up to the stellar surface. Figure 19 plots this as a function of the integration limit , for the ℓ = 1 modes used in Fig. 18. While some small-scale fluctuations resulting from the mode structure are present, the overall behaviour in all cases is a relatively steep increase in the most central regions, where the rate of growth is largely independent of . Further out where the configurations undergo their finer-scale spatial reversals, the integrals plateau to a constant value that is lower for larger . Since all curves are monotonically increasing (if one ignores the small-scale fluctuations), this means that the integrand is always positive and so we conclude that for larger , it must be that the integrand itself is smaller over more of the volume. With aid of the plots showing total field strength (Figs S4-S7), it is apparent that this is driven by the steeper drop-off in field strength with radial distance for configurations with larger . Whether this is a peculiarity of the Prendergast solution or likely to be a property of real stellar magnetic fields is unknown. However, it is to be noted that a separate physical argument for weaker overall field strength in regions possessing smaller-scale field structure is the higher rate of Ohmic dissipation.
The inertial-frame frequency shifts (including the Coriolis force) for all modes of Model C are shown in Fig. 20, for two values of . As might be expected, similar to the case for f , the g-dominated mixed modes are far more sensitive to changes in than the p-dominated ones.

Implications for asteroseismology
Multiplet asymmetries are a known consequence of the Lorentz force, contrasting the Coriolis force which by itself only produces symmetric multiplets in the limit of slow rotation. While it may be tempting to use this signature as an indicator for the presence or absence of core magnetism, there are several complications. Firstly, magnetism is not the only possible cause of asymmetry (cf. Deheuvels et al. 2017). Secondly, the results here show that it is possible for multiplets to appear symmetrically split, by a value close to that associated with pure rotation, even when the Lorentz force is comparable to the Coriolis force. This tends to occur for intermediate obliquities ( ∼ /4); see e.g. Fig. 12 and Fig. 15 (middle row). This result is important because it means that the detection of symmetric multiplets cannot alone be used to rule out the existence of a magnetic field.
A closely related concern is the fact that magnetic fields of nonzero obliquity are supposed to increase the multiplicity from 2ℓ+1 to (2ℓ +1) 2 peaks in the inertial frame. Given this, it might be similarly tempting to interpret the observed lack of additional multiplicity as evidence of the absence of a magnetic field. However, it is clear from  Fig. 9), under the default rotation and field strength listed in Table 1, for six different values of the obliquity angle. Short, thick, horizontal black lines joined by thin dotted lines illustrate the splitting process, from the unperturbed mode (leftmost line at 1 = 0) into a multiplet of 2ℓ + 1 modes in the corotating frame (middle group), and then further into (2ℓ + 1) 2 modes in the inertial frame (rightmost group). In the rightmost group, the various modes have been horizontally offset for clarity. The different greyscale shades correspond to the absolute value of the associated coefficient, with black being 1 (maximum amplitude) and white being 0. For comparison, horizontal coloured dashed lines mark the position of the non-magnetic rotationally split multiplet in the inertial frame.  e.g. Fig. 9 (bottom right) that most of the additional components can have quite low amplitudes, perhaps an order of magnitude smaller than the dominant components, even in situations where the Lorentz force is comparable to the Coriolis force and is substantial. Thus in practice, many of the additional peaks may not stand out well above the noise, and only the 2ℓ + 1 dominant components along the diagonal of the a matrix have sufficient amplitudes to be identified. It is therefore possible that magnetic fields of substantial strengths are present in stars for which this scenario has previously been ruled out on the basis of the lack of multiplet asymmetry and/or absence of enhanced multiplicity. Rather, more detailed modelling on a caseby-case basis may be required to support such conclusions.
There may be some specific behaviours which could provide useful arguments for or against the existence of a core magnetic field. As noted throughout Section 4, magnetic effects are more important for modes with shorter wavelengths and higher inertias. Signatures such as multiplet asymmetry and enhanced multiplicity are therefore expected to be most pronounced for modes that have lower 0 , higher ℓ, are more g-dominated, and in stars with greater central condensations. In particular, the 1/ 0 dependence of the magnetic shift and its tendency to systematically increase the centroid frequency of a multiplet is in contrast to the behaviour of rotation, which at leading order has no 0 dependence and does not alter the centroid frequency. Since g-dominated mixed modes are more heavily affected by a core field, the centroids of their associated multiplets would be displaced further from the unperturbed values compared to p-dominated mixed modes. For example, if one were to perform forward modelling of an evolved star to reproduce well the p-dominated mode frequencies while neglecting the existence of a core field, then this would systematically underpredict the frequencies of g-dominated modes. The discrepancy would be worse at lower 0 and higher ℓ. Demonstration of such behaviour or lack thereof could provide a convincing case for/against core magnetism, but may be tricky in practice owing to the need to detect sufficient numbers of g-dominated mixed modes.
As a comment, this property of the Lorentz force to produce systematically positive frequency shifts owes to the fact that the integral in (28) is positive definite as long as the integration limits extend beyond f (outside of which = 0). Upon integration by parts, the boundary terms may thence be neglected, and the integral can be approximated by ∫ |∇ × ( × B)| 2 d ≥ 0. Note that in this paper we choose to use the more general/exact expression in (28): this is necessary in particular for Fig. 19, since the radial integration limits are truncated within the field region.
Regarding seismic inference of the geometry/topology of a magnetic field, at first glance it may seem as though a near-infinite number of parameters would be needed to parametrise a field configuration of realistic complexity, making this a daunting prospect. However, this study shows how the problem may be made more tractable, since we have been able to parametrise a fairly complex magnetic field using just four global scalars: , A,cen , f and . Even without extensive characterisation of the g-dominated mixed modes, constraining some of these parameters may be possible using p-dominated mixed modes alone (cf. Fig. 16, for which the bottom row are p-dominated modes).

Limitations
The Prendergast solution may be a highly convenient modelling choice for the magnetic field, and supported by a large number of theoretical and empirical considerations, but ultimately it is difficult to know how accurately this would represent the field in the core of an actual star. As noted, for each stellar model there is actually a whole family of solutions, having different . There are suggestions Figure 15. Inertial-frame frequency shifts versus unperturbed frequency, for all ℓ = 1 modes of Models C (left) and D (right). Values of the field strength and obliquity are shown in the panel headers. Within each panel, the different colours correspond to the corotating-frame modes, which are each split further into 2ℓ + 1 modes in the inertial frame. The different symbols represent the different components: upward triangles represent = +1 (retrograde modes), circles represent = 0 and downward triangles represent = −1 (prograde modes). In addition, each symbol is filled with a colour whose saturation varies according to the value | |. Fully saturated (i.e. solid red/blue/black) corresponds to | | = 1, while white corresponds to | | = 0.
that the configuration corresponding to the smallest root is the one most likely to occur in reality, as it has the lowest associated energy (Duez & Mathis 2010). However, it may also be argued that the topology of the final configuration should depend on that of the dynamo field at the point of cessation, assuming helicity conservation during relaxation as would be expected in ideal MHD. Thus it may be considerably more complex than the smallest-configuration, although the substantial uncertainties associated with dynamo physics leave this issue open. In this work we have considered higher-configurations and found these to produce smaller magnetic shifts, for the same central strength. This appears to be a consequence of a steeper rate of fall-off in field strength with radial distance, which is in turn linked to a larger fraction of the volume having small/nearzero field amplitude when there are more spatial reversals. While this might not seem inconceivable in reality, and would be favoured by physical arguments involving Ohmic dissipation anyhow, it is to be borne in mind that ultimately these results are for a specific type of field configuration, which may have certain peculiarities.
It is noticed that the rotational frequency shifts considered for Model D are comparable to the mode spacings, so it may be quasi-degenerate perturbation theory that would more appropriate. However, work done by Loi & Papaloizou (2019) suggests that first-order theory remains reasonably accurate at predicting the frequency shifts in this regime; rather, it is the accuracy of constructing the perturbations to the eigenfunctions (via the coefficients ) that suffers the most. However, that is not of interest here, and so we adopt first-order theory throughout. It is also worth commenting that the chosen field strength for Model D, which was adjusted so that magnetic and rotational frequency shifts would be comparable, is near the estimated critical strength at which gravity-Alfvén wave resonance and mode conversion might occur. While the effects of this cannot be directly handled by perturbation theory, the work of Loi (2020a) found the existence of directions in which gravity wave propagation may not be significantly affected by the Lorentz force, and others where mode conversion ("trapping") occurs efficiently. This anisotropy predicts that sectoral modes (with respect to the magnetic axis) may be relatively unaffected, while zonal modes are preferentially damped, which would act to depress certain | | coefficients, possibly pushing them below detection limits and thus lowering the apparent multiplicity. Further consideration of this  possible "hybrid" approach for incorporating the effects of a dynamically significant field is beyond the scope of the current work, but may represent a possible avenue for future studies.

SUMMARY
In this paper we have presented a phenomenological study of the frequency splitting patterns induced by magnetism and rotation. We applied this to the case of evolved stars containing core magnetic fields, and examined how this is influenced by properties such as the obliquity and topology of the field. Realistic mixed configurations are used which are fully parametrisable by just four global scalars, demonstrating that seismic inference of the properties of a core magnetic field may be a reasonably tractable problem. It is found that greater radial complexity of the field (more fine-scale structure) leads to smaller magnetic contributions to the frequency shift, owing to a faster fall-off in field strength away from the centre. Topological variations in core fields among evolved stars, due to differing types of dynamos or otherwise, may therefore speculatively play a role in explaining some of the diversity in e.g. mode depression behaviour.
A common assumption is that the lack of multiplet asymmetries and/or lack of increased multiplicity implies the absence of a magnetic field. However, the results of this work show that it is possible for both of these to be exhibited even when the Lorentz force is comparable to the Coriolis force, i.e. magnetic and rotational contributions to the frequency splitting are comparable. This arises from the non-straightforward effects of obliquity and the nature of -mixing, underscoring the subtlety of the Lorentz force and suggesting the need for more complete modelling of rotational and magnetic effects to support such conclusions.
Distinctive signatures of magnetic splitting that may allow it to be disentangled from rotational splitting include its general 1/ 0 dependence, positive definite nature, and tendency to affect modes of small spatial scales the most (due to the -dependence of the  Table 1. Rotation and field strength are set to their default values, and = 0. Red, green and blue correspond to ℓ = 1, 2 and 3 respectively, while different symbols correspond to different values of ′ . In increasing order from ′ = −3 to +3, these are denoted by pentagrams, squares, downward triangles, circles, upward triangles, diamonds, and asterisks respectively. Alfvén frequency, in contrast to rotational effects which are independent of ). Hence these signatures would be more pronounced for higher ℓ, at greater central condensations, and among g-dominated mixed modes. It would be interesting to apply these concepts to analysing real objects in which previously, on the basis of fairly "normal"-looking splitting patterns, magnetic fields may have been inferred not to exist. Note that the tools used in this study are publicly available (see Data Availability section below); those interested are encouraged to appropriate these for their own purposes.

ACKNOWLEDGEMENTS
The author is supported by funding from Churchill College, Cambridge through a Junior Research Fellowship. We thank the code development team members for all their efforts put into their publicly available stellar evolution code (version r11701).

Source
codes and smaller data files can be found at the following link: https://github.com/STCLoi/topology-obliquity-paper/. For any queries, or access to larger data files, please contact the corresponding author.  Below are the components of J ′ ( , ′ , ′ ), the Eulerian perturbation to the current density, which appears in the expression (28) for the components of the matrix M mag . These are for the special case of the Prendergast solution (see Section 2.3 for details). Arguments have been suppressed for brevity, but recall that ℓ , ℓ and Ψ are sole functions of , and . Inertial-frame frequency shifts versus unperturbed frequencies for all modes of Model C, for two different values of , these being the second (top) and eighth (bottom) roots of (4). Colours and marker styles are used in an identical manner to Fig. 17.