Numerical modelling of imposed magnetohydrodynamic effects in hypersonic flows

Weakly ionised plasmas, formed in high enthalpy hypersonic flows, can be actively manipulated via imposed magnetic fields - a concept termed magnetohydrodynamic (MHD) flow control. Imposed MHD effects, within flows which exhibit multiple shock interactions, are consequential for emerging aerospace technologies: including aerodynamic control via magnetic actuation. However, numerical modelling of this flow type remains challenging due to the sensitivity of feature formation and the real gas modelling of weakly ionised, electrically conductive, air plasma. In this work, numerical simulation capabilities have been developed for the study of MHD affected, hypersonic flows, around 2D axisymmetric non-simple geometries. The validated numerical methodology, combined with an advanced 19 species equation of state for air plasma, achieves quantitative agreement between simulation and experiment for a Mach 5.6 double cone geometry with applied magnetic field. Numerical studies are conducted for varied conical surface angle and magnetic field configuration. This paper demonstrates how, for hypersonic flows with complex shock interactions, the MHD affected flow is not only augmented in terms of shock position, but may exhibit topological adaptations in the fundamental flow structure. A classification system is introduced for the emergent flow topologies. The applied numerical studies examine the mechanisms by which the magnetic field configuration influences the MHD augmented shock structure, leading to: (1) differences in magnitude of MHD enhancement effect, and (2) structural adaptations of the flow topology. Most critically, classes of conditions are identified which produce topological equivalence between the magnetic interaction effects and a generalised mechanical control surface.


I. INTRODUCTION
In high enthalpy hypersonic flows, the pressures and temperatures in the shock layer can become sufficiently high that dissociation and ionisation reactions produce a weakly ionised plasma. Owing to its elevated electrical conductivity, the plasma can be actively manipulated via electromagnetic effects when a strong magnetic field is generated from within the vehicle.
The MHD flow control concept is depicted in FIG. 1 for a generic blunt body. The charged particles travelling within the shock layer pass through the magnetic field lines (B), which induces a strong circumferential electric current (J). This electric current interacts with the magnetic field to produce a Lorentz force (J × B) which acts macroscopically on the fluid. The resultant effect is to force the bow shock further away (∆) from the vehicle surface (increasing the bow shock total surface area), and is termed shock layer enhancement 1 .
At present, prospective applications of MHD flow control for aerospace technologies are an active area of research. These include: the MHD by-pass engine (acts to accelerate ionised flow out of a combustor to generate additional thrust) 2,3 , hypersonic inlet control for scramjet engines (the leading oblique shock angle is magnetically controlled to achieve on-design specification at off-design flight condition) 4-6 and the prospect of active thermal protection via magnetic heat shielding [7][8][9] as well as drag enhancement 10,11 FIG. 1. MHD flow control concept: streamlines (V) are shown, carrying charged particles inside the bow shock. Two local tangential components of the larger dipole magnetic field are labelled B. The induced circumferential electric current J θ crosses through the magnetic field to produce a Lorentz force in the J × B direction, enhancing the shock stand-off distance by ∆.
for atmospheric reentry vehicles. In flows with insufficient ionisation, alternative methods are used to locally elevate electrical conductivity. Such methods (including plasm actua- tors in the boundary layer) often generate an additional electric field, which creates a supplementary mechanism in the Lorentz force 12 . In hypersonic flows where a sufficiently high level of electrical conductivity is generated via direct impact of the freestream with the vehicle body, the imposed magnetic field can directly affect the flow 13 -this is the focus of this work.
Extreme forces are experienced along the surfaces of hypersonic vehicles. Control flaps and surfaces require high power mechanical actuators and multiple moving parts to adjust surface angles dynamically during flight. MHD flow control presents the desirable possibility of replacing mechanical actuation of the control surface, however, such research is in the early stages 14 .
Experimental and numerical research for MHD control of hypersonic flows has generally focussed on simple geometries 15 which produce an elementary shock structure: such as the single bow shock formed around blunt bodies 16,17 and the leading oblique shock formed over wedges 18,19 and cones 15,20 . In the case of multiple angled surfaces -which realistically represents a mechanical control surface, such as a flap or ramp -the coupled fluid dynamic and magnetodynamic effects become complex. Shock wave boundary layer (SWBL) interactions lead to the formation of a separation region, and produces multiple shock interactions. With sufficient ionisation in the flow, the high conductivity low velocity regions which form within the steady state shock structure, are sensitive to manipulation via magnetic fields 21 .
The double cone is considered a useful test geometry since the resultant flow exhibits all of the relevant features of common 3D internal engine geometries or aerodynamic control surfaces. FIG. 2 depicts the comparison between shock struc-ture control via imposed MHD effects and traditional mechanical surface actuation.
In the case without magnetic fields, detailed experimental research has been conducted on the 2D axisymmetric double cone 22 . Even without electromagnetic interactions, high enthalpy laminar flows with SWBLI and real gas effects are challenging to numerically simulate, and have been the subject of many studies [23][24][25] .
In the case with a magnetic field generated from within the model, experiments were conducted by Wasai et al. 14 in an advanced expansion tube facility, for a hypersonic flow (with air as the test gas) over a double cone geometry. The experiments were able to demonstrate and measure an MHD enhancement effect due to the applied magnetic field on the flow which exhibits multiple shock interactions. Simulations of the experimental conditions were conducted by Nagata et al. 26,27 who were able to demonstrate similar qualitative effects as observed in the experiment.
A number of prevailing numerical modelling challenges have been identified for this flow type. These include: the balance between computational tractability and physical accuracy of the governing mathematical model (the coupled fluid-electromagnetic flow physics exhibits complex multiscale phenomena), realistic thermochemistry modelling for weakly ionised air plasmas, accurate resolution of electric transport properties (which is typically related to the thermochemistry model), and numerical solvers which are robust for the hypersonic regime over complex geometries, including the stable implementation of Lorentz forcing dynamics.
In the current work, 2D-axisymmetric simulations are conducted to investigate the physics of hypersonic flows over double cone geometries, with imposed MHD effects. Notably, the developed model includes a new extended 19 species equation of state (EoS) for air plasma applied to hypersonically formed plasmas in the equilibrium regime. This EoS enables the direct and accurate computation of electrical conductivity which is a key parameter driving the Lorentz forcing dynamics. The numerical solvers are implemented within the AMReX framework 28 , which permits block-structured parallelism and hierarchical adaptive mesh refinement, resulting in a high effective resolution at low cost and scalability of the computation. Efficient and robust embedded boundary geometry (via a rigid body Ghost Fluid Method) permits arbitrarily complex geometries to be simulated within the structured Cartesian mesh, and are shown to be stable for strong shockbody interactions.
The numerical model is validated against experimental results. The model is then applied to study the double cone configuration with imposed MHD effects. Previous MHD studies on simple geometries -blunt bodies, single cones and wedges -produce an elemental shock structure: bow shock or oblique shock. When the magnetic field is imposed upon the flow, this causes a self-similar MHD enhancement effect: congruent shock structure, but with a change of position. In these cases, the enhancement can be effectively quantified through the measured change in one variable: % increase in shock stand-off distance along the stagnation line for a blunt body, or the increase in angle of an oblique shock wave for conical or wedge geometries. The numerical studies of this work, however, demonstrate how, for a complex flow field combined with an imposed magnetic field, the base flow and the MHD affected flow can become structurally dimorphic. That is, the topology of the magnetically augmented flow is often not selfsimilar to its respective base flow case.
Formational differences in the flow structure are critical within a context of aerodynamic prediction for flight control. Therefore, this work examines the flow physics of the complex, and coupled, magneto-fluid behaviours with a focus on the conditions and underlying mechanisms which lead to structural variation of the steady state flow fields. In terms of the magnetic actuation control surface concept, not only the measurement of MHD shock enhancement effect, but the identification of structurally similar flow morphologies between MHD actuation and mechanical surface actuation has technological significance.

II. MATHEMATICAL MODEL
The magnetohydrodynamic behaviour in a compressible flow regime exhibits interdependent phenomena across disparate scales. For numerical solution to be tractable, assumptions and simplifications are necessary.
In the case of aerospace plasmas -where collision frequency within the plasma is sufficiently high, flow length scales are significantly larger than the Debye length, and time scales are larger than the reciprocal of the plasma resonant frequency -a continuum magnetohydrodynamic (MHD) model governs 29 .
Broadly across the literature, hypersonic plasmas are con-sidered in the low magnetic Reynolds number regime 13,30,31 where Re m is of the order of unity or less. Under this regime the imposed and induced magnetic fields are sufficiently decoupled for the static imposed magnetic field to be unaffected by the flow.
Where σ 0 is the nominal fluid electrical conductivity, V 0 the total velocity, L 0 the characteristic length, ε 0 is the air permittivity, and c the speed of light in a vacuum.
The governing system equations are derived by combining the Navier-Stokes equations for fluid dynamics (accounting for fluid viscosity and thermal conductivity) with Maxwell's equations for electromagnetism. The derivation under the low Re m weakly ionised plasma regime results in resistive MHD equations which account for Lorentz forcing and Joule heating.
An important dimensionless number throughout this work is the magnetic interaction parameter (Q MHD ) defined as: Where B 0 is the peak magnetic field strength. The ability of the imposed magnetically field to dynamically affect the flow is best determined by Q MHD .

A. Governing system equations
For studies of axisymmetric vehicle geometries, the axis of symmetry can be exploited, converting the 3D problem into a 2D plane in cylindrical coordinates. Therefore the MHD system equations are given in r-z space as: U, F and G denote the state vector, and r-z directional inviscid flux vectors respectively, and S c denotes the associated cylindrical source term vector. Similarly, F v and G v denote the r-z directional viscous flux vectors and associated cylindrical source term S c v . S MHD denotes the MHD source term vector.
For the 2D r-z system, we denote the velocity vector in the (r,θ ,ẑ) direction to have components (u, w, v), and so the system vectors are defined as: where total energy E is given by: E = ρe + 1 2 ρu 2 . Assuming a Newtonian fluid, the viscous stress tensor is given by: where µ is the fluid dynamic viscosity and using stokes hypothesis λ = − 2 3 µ. ζ is the fluid's thermal conductivity. The flow is assumed to be laminar, and dynamic viscosity is computed via Sutherland's law. Thermal conductivity is computed with constant Prandtl number = 0.71 for the ideal gas case, and ζ is obtained from the tabulated (precomputed) state data of the EoS plasma19X (described in the next section).
Then the viscous flux vectors are derived as follows: u r ) µw z uτ zr + vτ zz + wτ zθ + ζ T z      and the following viscous source term vector results from the cylindrical system conversion: Some additional assumptions and resultant simplifications can be made in defining the S MHD vector.
The full expansion of Ohm's law includes a Hall current contribution. This additional current drives the generation of an electric field which has a complex interaction with the gas dynamics -termed Hall effect. However, the Hall current can be neglected when the Hall coefficient is sufficiently low. Additionally, several authors have shown that for a conductive body, the Hall effect diminishes the MHD flow control, whereas, in the case of a perfectly insulated body, the overall flow structure is essentially unaffected by the Hall effect 8,32,33 . Regarding the leading surfaces as an insulated material (via the electrodynamic boundary condition at the body), studies of MHD flow control justifiably neglect Hall effect 7,34? , therefore applying the generalised Ohm's law as given by equation 5: The electric field E is a conservative field where E = −∇φ and φ is the electric potential. Therefore the governing system is comprised of a system of 5 mixed hyperbolic-parabolic type PDE's as defined within equation 3, coupled with an additional Poisson-type elliptic PDE maintaining divergence-free current density: The magnetic field is a static imposed dipole field with dimensional components B = (B r , B θ , B z ), and defined analytically via: where r is the radial vector from the dipole centre, and m is an orientation vector aligned parallel to the dipole centreline. Within the axisymmetric governing system where r and m are defined in the r-z plane, the B θ magnetic field component is zero.
For electrodynamic boundary conditions which treat the surface as insulated, equation 6 can be solved analytically as φ = 0, which means the elliptic current density equation does not require numerical solution and electric current density is given by: Therefore the Lorentz forcing terms in the MHD source term vector can be defined directly, rendering the computation as vastly more expedient:

B. Equation of state for air-plasma
The system equations are closed with an advanced 19 species equilibrium air-plasma equation of state (EoS), applied for the first time to applications in hypersonically generated air plasmas. The new 19 species approach to the thermochemistry modelling is an extension to the standard 11species model 35 , and includes a direct and accurate computation the electrical conductivity transport property, which has previously posed challenges in MHD studies 36 .
The advanced EoS considers 19 species reaching chemical and thermal equilibrium, and therefore can be applied to aerospace plasma problems determined to be in local thermodynamic equilibrium (LTE). The EoS is an extension of the precomputed plasma19 database originally developed by Träuble in 2018 37 (which utilised key theory from D'Angola 38 and Colonna 39 )), and as published in 2021 for the application of simulating air-plasma generated during lightning strike 40 . Extensions were made to this database in this work, based on the original theory but in order to generate properties over a wider density-pressure range necessary for the simulation of air-plasma generated by hypersonic flows at the low pressure and density freestream conditions typical of this regime. Additionally, an extension to the computation of electron number density was implemented, following the theory of D'Angola 41 . The new extended EoS will therefore be referred to as plasma19X.
The 19 species considered in the air-plasma are: The development of the 19 species EoS considers a large number of basis chemical reactions, where every such reaction is governed by a law of mass action corresponding to dynamic equilibrium, and total mass is conserved across all species. For a specified intensive property pair (e.g. pressure, temperature) a large, determined system of non-linear equations governs the precise species composition. Numerical solution for species composition, and resultant thermodynamic properties from detailed curve fits which follows the work by D'Angola 38 . For the remainder of the specifics of the developed theory for plasma19, please refer to the detailed description in the documented work 37 .
Following the theory of plasma19's development, thermodynamic properties are physically extended over the low pressure and temperature range required for freestream conditions of hypersonic flight cases. A method to very accurately compute electron number density down to very low air plasma temperatures was implemented for the extended tabulated plasma19X EoS following the extended theory proposed by D'Angola 41 . This is critical to the plasma properties of the temperature band of interest in hypersonically formed plasmas (approximately 2,000-10,000 K). The extended plasma19X EoS can compute thermodynamic properties, transport properties and species composition for a density range of 10 −5 kg/m 3 -10 kg/m 3 and a pressure range of 7 Pa to 18.23 MPa, which corresponds to a temperature band of approximately 10 K -60,000 K. The state properties computed by plasma19X have been extensively validated against empirical measurements 40 .

C. Electrical conductivity model
Accurate modelling of any magnetogasdynamic effect is underpinned by the accuracy of the computed electrical con-ductivity. This is because conductivity drives the Lorentz forcing dynamics of the flow: As Bisek 36 highlights, accurate conductivity modelling typically involves semi-empirical approximations and remains problematic in MHD studies.
Ideal gas models (widely used in MHD studies 8 due to the computational efficiency) cannot inherently compute electrical conductivity, and so, adopt a highly simplified supplementary model. The supplementary conductivity model is typically an analytic power function of temperature 30,42,43 : where constants σ 0 and T 0 are manually fitted to the problem, and exponent n takes a value in the range: 0-4 depending on the test case.
Real gas models (both equilibrium and non-equilibrium) also adopt supplementary conductivity models which are functions of electron number density (N e ), and therefore critically depend on N e as computed by species dissociation and ionisation of the full thermochemistry model. A standard form for weakly ionised equilibrium plasmas is given by: where n e is the electron number density, e is the electron charge, m e is the electron mass, and ν coll is the effective electron-neutral and electron-ion collision frequency. Several different formulations for ν coll (where the gas is assumed to follow the Maxwell Boltzmann distribution) are commonly considered 50 . These are typically in the form of analytic functions of pressure and or temperature, where constants vary depending on the empirically identified regime 36 .
The electrical conductivity given by plasma19X, however, was computed by the detailed theory of perturbative Chapman-Enskog method, which derives a full set of collision integrals from the species compositions, and a third order approximated diffusion coefficient (D ee (3)). which then computes electrical conductivity as: where n is the total number density and k B is the Boltzmann constant. The detailed derivation can be found in Capitelli et al. 51 .
FIG. 3 demonstrates the accuracy of the computed σ by plasma19X as compared with experiment and other commonly employed supplementary models. As can be seen in FIG. 3, the non-physical species molar fractions over the low temperature range, and resultant anomalous electrical conductivity, render EoS11 as problematic for LTE MHD studies. All  45 and Yos (1963) 46 where a low-temperature equilibrium plasma was generated via a stabilised DC arc applied to a column of air maintained at atmospheric pressure. Otsu 8 and Bush 47 use the simple analytic power law approximation defined in equation 10 with exponents 2 and 4 respectively. Solution to the Boltzmann equation is as per the solver developed by Weng and Kushner 48 for two different levels of resolution. EoS11 is an 11-species equilibrium model as implemented by Villa et al. 49 . models must be both smooth and accurate to produce stable and realistic computational results for MHD studies. While different models perform well in different temperature bands, this is often achieved by manually fixing parameters to the computed temperature and profile of the flow. Plasma19X is found to be smooth and accurate (without any manual input), especially over the temperature band of interest for hypersonically formed plasmas in the equilibrium regime.
The accuracy of plasma19X can be attributed to the detailed method of derivation of σ in equation 12 including the additional species dissociations and ionisations captured at relatively low air plasma temperatures in the 19 vs 11 species approach. With relevance to the application of this work, the double and triple ionised species captured by plasma19X are only present at higher temperatures (>20,000 K), and largely do not effect the state in the temperature band of interest (approximately 2,000 K -10,000 K). The negatively charged O − , O − 2 and N − 2 species however, are not captured by an 11-species model but are computed to be present in small concentrations in the low air plasma temperature range. The implementation of electron molar fraction is accurate down to temperatures as low as 50 K as compared with verified result of D'Angola 41 . Plasma19X also has the advantage that σ is obtained directly via the pre-computed tabulated EoS, rendering the model as computationally efficient. The numerical model of this work is developed within the AMReX framework 52 . AMReX provides a platform for hierarchical adaptive mesh refinement (AMR) and massively parallel computation of domain sub-blocks. All numerical methods were implemented directly within this framework in developing the model of this work. The hierarchical AMR (depicted in FIG. 4) and subdivision of the regular Cartesian mesh, which covers the full domain, results in efficient and highly scalable computation.
Due to the full domain discretisation of the AMReX framework, geometry is implemented via an embedded boundary method. A sharp interface rigid body Ghost Fluid Method (GFM) was implemented, such as that employed by Bennett et al. 53 , based on the approach of Sambasivan and UdayKarmar 54 .
By implementing embedded boundary geometry via a rigid body GFM, the mesh is not subject to distortions and skewing which sometimes occur in the case of body fitted meshes, especially as geometries become more complex. The Riemannbased GFM is robust for very strong shocks which interact with the embedded surface. The adaptive refinement approach ensures high resolution is achieved automatically and locally in the vicinity of high density gradients and at the body-fluid boundary. This combination of methods therefore permits the simulation of non-simple geometries, in the presence of strong interacting shocks, with multiple levels of hierarchical AMR applied to achieve high effective resolution and solved efficiently via MPI-based computation with high scalability across many nodes.

B. Numerical solver methodology
The combination of algorithms, implemented within the AMReX framework, is a recently developed methodology, outlined in a previous publication 55 , designed for the efficient simulation of hypersonic MHD systems.
The inviscid system fluxes are evaluated numerically via a complete wave approximate Riemann solver (HLLC) with high order extension (MUSCL) which is second order accurate in space and time. Viscous derivative terms are evaluated using the central differencing method at cell interfaces with divergence then computed through the finite volume construction. The transient solver can be used for time accurate simulation of transient test cases, and also reaches a stable and convergent steady state, for such test cases.
The system is solved via an operator split strategy whereby the hyperbolic and diffusive fluxes are computed successively with the result of each time-wise integration input as the initial condition for the subsequent operation. Fluxes are also dimensionally split in this manner. Cylindrical source terms are integrated via Strang splitting to maintain second order accuracy.
The Lorentz forcing ODEs of the MHD model can be expressed in the condensed form: where η = σ −1 is the plasma resistivity. An explicit Euler integration approach is adopted, where an additional momentum half-step is utilised: The momentum half-step, and resultant u n+ 1 2 are carefully constructed to ensure the J × B specifically effect momentum whilst conserving internal energy. The joule heating term η n (J n ) 2 only acts upon system internal (thermal) energy. See full derivation in Villa et al. 49 .
Some careful treatment of system stability is required for dimensionally and operator split systems when computed within an AMR framework.
The stable time-step is restricted as the minimum of the hyperbolic and diffusive stable time steps, given by: where ∆x d,i is each dimensional spatial step size. ∆t di f f is based on the formulation recommended by Gokhale 56 . Under the operator split solution strategy, the explicit time-marching evolution is stable in 2D (including all source terms) for CFL=0.9. However, stable time-step is computed on the 0 th (coarsest) level of the Cartesian mesh, with higher levels subcycled with matched integer fractions of the base-level timestep. Therefore, in consideration of multiple AMR levels and the combined hyperbolic-parabolic nature of the Navier-Stokes system equations, the stable time-step as computed at the 0 th mesh level is: where N_MAX is the highest AMR level and a refinement factor of 2 is assumed.

IV. VALIDATION TESTS
A large number of tests have been performed to extensively validate the new numerical model. A set which collectively validates the key physics of this study is presented here.

A. Low Re m MHD system
The accepted benchmark validation test for a hypersonic MHD system under low magnetic Reynolds number assumption, is the theoretical test originally proposed by Poggie & Gaitonde 42 . It has been replicated subsequently by Damevin & Hoffmann 57 and most recently by Fujino et. al 58 . For the purpose of validating the underlying resistive MHD model, the model of this work is reduced to its inviscid Euler equivalent, closed with an ideal gas equation of state for consistency with this validation test cases in the literature. For the simplified ideal gas model, the electrical conductivity is given by the analytic power law of equation 10, where σ 0 = 300 (Ωm) −1 and n=0 which therefore sets a constant post-shock conductivity through the shock layer (initiated via the condition: T > T ∞ ). B 0 is the maximum absolute magnetic field strength outside of the vehicle, and the characteristic length R 0 is the radius of the spherical body.
The MHD enhancement effect is shown to increase shock stand-off distance with increased Q MHD as defined in Equation 2. For this test Q MHD is varied from 0-6 by corresponding maximum B-field strength. 7 is typically chosen due to the availability of smooth experimental data, and results of multiple numerical studies for comparison. Since the test gas is non-reactive nitrogen (N 2 ), the ideal gas EoS was used to close the system. Very good agreement between the simulation, experiment, and the results of other numerical models is demonstrated in Figure 9. To show the effect of increasing the magnetic field strength (as per the Nagata simulations), a simulation is run with a stronger B-field: fixed at the same dipole centre location but scaled to a larger magnitude corresponding to 1.5 T at the reference point.
The MHD effects are depicted by the overlaid pressure contours for different magnetic field strength in FIG. 11. The small but clearly observable magnetic augmentation of the shock position in the B = 0.36 T case demonstrates the importance of high resolution numerical methods to sharply capture the shock wave position.
The experimental work resolves a property trace through lines parallel to the z-axis of symmetry at different radial distances (r = 4mm, 8mm and 10mm) along the model to determine shock position more precisely. These brightness intensity property traces along the different radial lines are shown in FIG. 12, and the increase in leading shock position due to MHD effects of the imposed magnetic field can be measured from these graphs. The leading shock position is considered to be where the property trace of brightness intensity first rises in value. Since multiple runs are presented in the one graph, the minimum offset and maximum off-set between the no-B and B = 0.36 T cases of all test runs is measured and a % increase in the shock stand-off distance is presented as the range of values in Table I.
The important conclusion from the experimental results is that there is a clearly observable, and measurable, MHD enhancement of shock-stand-off distance. From the pressure traces along the radial lines, this observable increase in stand- off distance is presented. In Table I, the measured increase in leading shock stand-off distance from the simulations of this work is compared to the measured increase from the trace taken along the radial lines in the experimental results published by Wasai et al. 14 . The results show an increase in shock stand-off distance which falls within the range measured from the experimental results. Therefore the predictions of this work are quantitatively in-line with experimental measurements for matched magnetic field strength. The agreement in shock stand-off enhancement indicates that the model of this work is able to realistically capture the complex MHD flow control effect. This is a key result for the complex flow field: a hypersonic flow with SWBL interaction over the double cone model coupled with magnetic interaction.

D. MHD effects verification
Results of the Nagata et al. numerical simulations showed that for a magnetic field strength matched with the experimental conditions, the MHD effects were negligible in the computed steady state flow field solution (see property traces of FIG. 12). However, the experiment does in fact reveal a small but observable increase in shock stand-off distance due to MHD flow control. By initialising the magnetic field strength in the numerical simulation to be much higher than that of the experiment, Nagata et al. were able to demonstrate and explore the qualitative effects an imposed magnetic field has on the flow, which agrees with the qualitative effect observed in the experiment.
As such, in this verification section we seek to compare the qualitative effects observed for the introduced strong magnetic field case (B re f = 1.5 T) with the effects reported in the previ-  14 . Middle) pressure traces taken at the same radial lines from the numerical result of Nagata et al. 26 . Bottom) pressure traces taken at the same radial lines from the numerical result of this work. Experimental and Nagata traces are presented as total distance along the central z-axis (labelled X-axis in the experimental results), and simulation traces are shown more clearly as distance along the z-axis measured directly from the surface intersection point as z=0 for each case. ous numerical study.
A comparison of computed electrical conductivity is also shown for the two cases in FIG. 13: base flow without magnetic field, and the strong magnetic field case. This plot more clearly shows the qualitative differences between the base flow and an augmented flow under strong magnetic influence. The strong B-field strength result agrees in terms of qualitative phenomenology with the simulation results of previous studies 26,27 . Specifically, we observe the MHD effects cause the shock layer over the second plane of the model to lift away from the surface, and for the stronger B-field the shock wave begins to become bow shaped. The separation point moves upstream, and shock triple point becomes extended. The electrical conductivity is a dominant parameter driving the computed Lorentz forcing dynamics, and FIG. 13 shows the enlarged extent of high conductivity plasma in the detached shock region of the B re f = 1.5 T test case.
The Nagata studies were successful in characterising the qualitative effects predicted by strong imposed magnetic fields on the flow augmentation. They concluded in their work that the two key factors limiting the accuracy of the model were: effectively capturing the sensitive feature formation at sufficiently high resolution where shock waves interact, and the accuracy of the thermochemical model (which underpins prediction of electrical conductivity). As such, the model of this work aims to build on these identified factors by utilising high resolution, high efficiency, numerical methods combined with advanced thermochemistry modelling for air plasma and its electrical properties.

V. NUMERICAL STUDIES OF HYPERSONIC FLOWS WITH IMPOSED MHD EFFECTS
The validated numerical model is now applied to varied configurations of the hypersonic double cone geometry with an imposed dipole magnetic field.
Within the motivating context of hypersonic flight control, even small changes in shock position can affect the flight dynamics. For a non-simple shock structure, morphological adaptations -changes of the actual topology of the shock structure -are particularly important to identify and to classify. In relation to the surface actuation concept, the identification of topological changes in shock structure determines whether MHD actuation can morphologically replicate mechanical surface actuation effects.
These studies seek to explore and explain the effects of magnetic field configurations, combined with variations in geometry, in terms of the quantitative (% enhancement of detached shock position) and structural adaptations (different shock structure topologies) observed in the emergent flow fields. Through such studies, the conditions and mechanisms of action in the coupled fluid-dynamic and magneto-dynamic flow physics seek to be better understood. Conditions which amplify the MHD enhancement effect, and conditions which produce flows of morphological equivalence between surface geometry and magnetic actuation can be identified.

A. Shock structure classification system
The earliest experimental and theoretical work on the classification of shock interactions was conducted by Edney who examined the mechanisms of interaction between an oblique shock wave and a detached bow shock around a cylinder 61 . The emergent shock interactions were classified as Types I- VI. This work was then extended by Olejniczak et al. to determine inviscid shock interaction classifications when applied to double wedge geometries 62 . Though the shock structures differ between the cylinder and double wedge cases, corresponding shock interactions for the steady shock structures which form for an inviscid 2D flow over a double wedge can be determined for Types IV, V and VI. Olejniczak et al. also identify an additional shock interaction Type IVr which occurs over a double wedge but which does not occur in the original cylindrical shock classification system.
The corresponding Type VI and Type V shock structure classification of the local interaction near the kink for an inviscid flow over a double wedge is depicted in FIG. 14. As can be seen, when neglecting viscous effects the shock interactions involve an attached shock from the kink point of the double wedge geometry. However, this entire region becomes separated in the case of viscous hypersonic flows, therefore forming entirely different types of shock interactions.
Olejniczak et al. note in a subsequent study that the problem becomes more complex when viscous effects are accounted for 63 . The viscous SWBL interaction which leads to the formation of a separation region, produces a separation shock which alters the basic structure of the flow as it interacts with the other flow features. Additionally, growth of the separation region can cause the flow to become unsteady under certain freestream and geometric conditions 64 . Given the divergent phenomenologies, Olejniczak et al. conclude that it is not always possible, nor appropriate, to relate viscous double wedge flows to Edney's classification system. More recent work examines the complexities of viscous double wedges and control surfaces 65,66 , however, flows in this regime do not adhere to existing shock interaction classification systems, and no additional shock classification system is introduced.
This work considers viscous axisymmetric hypersonic flows over double cone geometries, as distinct from the double wedge studies. The dimensionality of the cone vs wedge promotes steadiness of the flow. Classification of the emergent steady state flow types serves as an important basis for analysis. Therefore, 3 classifications are proposed for the emergent flow topologies of this work. As distinct from Edney's classification system, and the double wedge equivalences, the classification system is not based on the singular shock interaction arising near the kink point, but the form of the overall shock structure as determined by the presence, absence or duplication of salient flow features. Specifically: the formation of shock triple points, and bounded subsonic regions. Therefore the topologies are delineated as classes, which consider all shock interactions (multiple may occur) within the total shock structure.
This classification system is shown in FIG. 15 and will be used in the analysis of this section. CLASS A represents the topology produced by hypersonic double cone tests in the CUBRC validation tests, where a single shock triple point, detached shock, and bounded subsonic region forms. CLASS B identifies flow structures where features surrounding the detached shock region ahead of the second surface become duplicated. Here a second shock triple point is associated with the formation of a second discrete subsonic region and a distinct secondary portion of the detached shock. CLASS C collectively identifies cases where the reflected flow over the second surface produces a shock inflection point, rather than a shock triple point, and no isolated subsonic region forms beyond the boundary layer and separation region. The defini- To first examine the effect of varying geometry on the base flow (no magnetic field), we conduct a parameter study of varied second surface angle (θ 2 ) . Inclination of the second conical surface generically replicates the effect of a mechanical control surface. It is identified from the previous studies that magnetic interaction is dominant in regions where the flow is at high temperature (high conductivity) and low velocity (low inertial forcing). Lorentz forcing has a linear dependency on |u| (see equation 19) whereas fluid inertial forces scale with |u| 2 . Therefore the magnetic interaction parameter: ratio of Lorentz forcing to fluid inertial forces is Q MHD ∝ |u| −1 . Increased flow stagnation and elevated electrical conductivity occur concurrently by increasing the angle of the second conical surface. We therefore conduct a parameter study with θ 2 surface angles of: We now introduce the dipole magnetic field in order to study the imposed MHD effects on the flow cases. For consistency with the experimental configuration of Section IV C, the magnetic field is configured to match the strong B-field case, with the dipole centre located at 19.5 mm from the conical point. The corresponding peak |B| value in the fluid domain is 2.25 T, occurring partially along the second surface, as shown in FIG. 16 (B 0 of equation 7 to produce this field is cited in figure).
Varying the second surface angle affects both the maximum |B| in the flow field and the angle of the B-field lines relative to the surface (and flow field). There is no perfect solution to isolate the effect of changed surface angle within the dipole magnetic field, given this is a complexly coupled problem. However, configurations can be achieved which closely control B-field orientation and magnitude relative to the actuation surface. The selected approach is to maintain the dipole centre at Z dist = 19.5 mm from the conical point, whilst fixing the shoulder point and kink point of the geometry on the r-axis, and extending along the z-axis as the θ 2 angle changes. This produces a magnetic field which is approximately congruent for each test case angle and B 0 is scaled to maintain precisely constant B max of 2.25 T, in the flow field.
FIG . 18 shows the steady state pressure contours for the base flow case as well as the enhanced shock structure of the MHD affected case for each conical surface angle. Reference lines are added to aid the analysis of MHD enhancement effect. Two classes of effects are observed as a result of the imposed magnetic field: 1. Alteration of the flow structure classification 2. Quantification of the magnitude of shock enhancement effect between cases Notably, for the θ 2 = 65 o case, addition of the magnetic field is observed to alter the steady state flow topology from CLASS B to CLASS A. The magnetic enhancement of the separation region alters the fluid dynamics downstream along the second surface, while Lorentz forcing acts to bow the shape of the detached shock in all cases. Combined, the resultant steady state forms to become CLASS A under the imposed MHD effects.
With reference to the guide line running tangential to the second surface (off-set at 2mm), for increasing θ 2 values: all 3 inclined surface cases become more bow shaped and show a decreasing segmented arc which exceeds the tangential guide line. Measurement of shock stand-off and % enhancement from the base case is calculated by taking the peak stand-off distance perpendicular to the second surface tangential line for both the base case and MHD affected case. Values are compared in Table II. The base flow conditions depicted in FIG. 17 indicate increased θ 2 promote conditions which would amplify the shock enhancement effect once a magnetic field is introduced. The counter-intuitive result which arises, is that the peak standoff distance and % enhancement both diminish for increasing θ 2 . The MHD enhanced shock structure lifts maximally away, perpendicularly from the surface lines, for the case of minimum θ 2 = 55 o . The separation region, however, becomes significantly more enhanced with increasing surface angle. In the base case flow (without magnetic field activation), the separation region is larger, the more inclined the second surface. The subsequent expansion of the separation region due to magnetic interaction results in maximum enhancement for largest θ 2 .
A plot of magnetic interaction parameter (scaled logarithmically from 0.1-10) is compares the different conical surface angles in Figure 19. The significant expansion of the separation region can be explained by the strong magnetic interaction which is observed to occur in the high conductivity low velocity separation region near the kink. Expansion of the separation region causes the separation shock to become more detached and inclined relative to the first surface. This causes the jet stream to rise and reflect further along the second surface (peak surface pressure point). Enhancement of the separation region diminishes the volume of fluid within the detached bow shock region where peak Lorentz forcing and magnetic interaction occurs in this section of the flow field. The peak Lorentz forcing vectors shift from the subsonic detached shock region and towards the separation region as the surface angle increases. This coupled effect between the detached shock layer and the separation region leads to smaller maximum stand-off distances in the upper bow shock region for increased θ 2 .  We now examine the effect of changing the magnetic field configuration on the MHD affected flow.
Moving the magnetic centre position along the z-axis maintains magnetic axisymmetry whilst altering the orientation of the magnetic field lines in the region where they emanate beyond the vehicle surface. Simply moving the dipole centre location with fixed B 0 in the dipole equation means the maximum B-field strength interacting with the fluid also varies. Note that |B| diminishes cubically with radial distance from the dipole centre. To best isolate the effect of B-field ori- Most critically, we observe that the Z dist = 30.0 mm case leads to the emergence of the CLASS B flow structure (where the non-magnetic base flow was of CLASS A). Further, for the case where the magnetic field is at position Z dist = 10.0 mm, we observe the emergence of the CLASS C flow structure. The schematic for the CLASS C flow structure is given in FIG. 15. As per the classification schematics, the flow is characterised by a shock inflection point, rather than the formation of a triple point, ahead of where the flow impinges on the second surface. This results in a fully supersonic postshock layer within the detached shock region.
Examining the general shape of the shock structures of FIG 21, the leading shock front over the second surface lifts away from the vehicle becoming broader above the shoulder for increasing Z dist , until the Z dist = 30.0 mm case, where we observe the emergence of the CLASS B flow structure. Since the shock wave splits midway along the second surface to form a second shock triple point, in the region centred at approximately the 6 mm radial line, the distance of the leading shock from the vehicle surface does not behave predictably with magnetic dipole position. Of the tested B-field strength and dipole positions, the Z dist = 19.5 mm case creates the largest shock enhancement effect over the region directly forward (z- Also, notably, the emergent CLASS C flow structure of the Z dist = 10.0 mm case results in a more suppressed shock structure. Compared to the base flow case, the detached shock remains closer to the second conical surface. Therefore, this magnetic field configuration could be considered to exhibit an MHD diminishing effect. As shown in the Q MHD with Lorentz forcing vectors plot, the peak forcing and magnetic interaction both occur close to the kink point. MHD forced deceleration of the flow at the kink causes the separation region to enlarge (see leading pressure contours of the separation region have migrated upstream along the first surface) which causes the flow to be deflected away from the kink at an earlier separation point, resulting in a separated shock which becomes reflected by the second surface further downstream than the base case. Since the flow behind the separated shock is more tangentially aligned with the incline of the second surface (as compared with the base flow), the shock is weakly reflected without forming a shock triple point. Absence of the triple point also eliminates the small stagnation region which forms just downstream of this shock feature (stagnation region retained in all the other flow cases). Note that the colour-mapped streamlines reveal a much higher (and supersonic) velocity magnitude behind the shock reflection point than the shock triple point of the other test cases. Greater velocity magnitude in this region diminishes Q MHD via two concurrent mechanisms: the non-stagnated gas results in a lower static temperature and FIG. 22. Examining the effect of magnetic field location on separation region: steady state pressure contours in black with computed velocity streamlines (coloured by total velocity magnitude) shown for each magnetic dipole centre location. thereby lower resultant electrical conductivity, while simultaneously the fluid maintains greater inertial force. These aspects dually act to diminish magnetic interaction. Recall that: Q MHD = (σ (p, T )|B| 2 R 0 )/(ρ|u|), and the flow exhibits both diminished σ (T ) and increased |u|. Consequently, as there is almost no magnetic interaction downstream of this reflection point, there is no MHD forced deceleration to bow the shock wave, and the shock position retains a flattened profile along the second surface. Therefore, counter-intuitively, the imposed B-field actually serves to diminish the overall width of the shock layer over the second surface.
To summarise, the CLASS C structure forms primarily since the peak B-field acts close to the kink point, enhancing the separation region without causing additional flow deceleration within the detached shock region. Combined, this results in a weak shock reflection instead of a shock triple point which flattens the downstream shock structure. In the Z dist = 13.5 mm test cases, the peak magnetic field value is still located very close to the kink, however, in the previous Z dist = 10.0 mm case the magnetic field lines were almost perpendicular to the z-axis at the kink, whereas in this dipole position the magnetic field lines emanate at approximately 40 o away from the z-axis at the kink point, and the magnetic field diminishes more gradually at this radial distance from the dipole centre (refer to FIG. 20). As a result, the Lorentz forcing vectors are also rotated. They now align to be almost direct deceleration vectors near the separation region (vectors run parallel to streamlines, forcing in the opposing direction to the fluid velocity), and the separation region becomes more enhanced in this test case. FIG. 22 shows the vortical streamlines of the separation region are notably enhanced from the base case flow (high Q MHD can be seen in this region of high conductivity and low momentum) which in turn makes the separation shock more detached and lifts the flow reflection point further along the inclined second surface. For this magnetic field position significantly greater Lorentz forcing occurs along the second surface. Lorentz forcing augments the flow further outwards from the second surface resulting in the formation of a shock triple point and subsonic downstream region which further enhances magnetic interaction due to the low flow momentum and high conductivity. As a result we see a notably bowed shock structure lifted away from the second surface, predominantly in the region just downstream of the shock triple point. The position of the shock wave above the shoulder eventually converges closely with the Z dist = 10.0 mm case, where magnetic interaction is weak or absent.
The Z dist = 19.5 mm case exhibits the same primary dy- namics. However, with peak |B| shifted further down the second surface, the magnetic interaction and Lorentz forcing is more significant in the detached shock region. This results in a more bowed shape of the detached shock and prolonged stem of the shock triple point. This creates the coupled effect of increasing the size of the subsonic high-conductivity region, which amplifies the total shock enhancement effect within the detached shock layer. CLASS B flow structure: Z dist = 30.0 mm case By placing the magnetic dipole centre far from the leading conical point, the peak |B| emerges just above the vehicle shoulder, and the magnetic field diminishes more gradually in intensity. The region of the dipole magnetic field that interacts with the flow, exhibits field lines which are more homogeneously oriented (less variation in the angles of the B-field vectors across the domain). The resultant flow field is very interesting for this test case. The Lorentz forcing acts to directly decelerate the flow velocity over the vast majority of the second surface shock layer (forcing vectors shown to almost directly oppose flow motion). The strong MHD deceleration forces which persist in the region around and above the shoulder effectively creates torque on the leading shock structure, which accumulates in angular force as the shock stand-off distance increases above the shoulder (torque is larger for forcing at a greater displacement), ultimately leading to a split in the detached shock and a second shock triple point being generated. This significantly alters the resultant shock structure (formation of CLASS B from CLASS A base flow), as seen through the comparison of pressure contours (FIG. 21).
A Mach number plot of the Z dist = 19.5 mm and Z dist = 30.0 mm cases are compared in FIG. 24. The split shock structure has a number of consequences: a second supersonic jet stream forms in the second surface shock layer, as well as a large second subsonic region forming downstream of the second shock triple point. The consequence of the shock splitting is that MHD enhancement effect is diminished (in terms of % increase in shock stand-off distance) along the first portion of the second surface, but remains significant above the shoulder. Precise prediction of magnetic dipole position and B-field strength combinations which reach the angular force threshold for shock splitting would be of high importance for design in MHD flow control.

All cases
In terms of maximising the MHD enhancement effect for a given magnetic field strength, the orientation of the magnetic field lines are best optimised when oriented perpendicularly to the fluid flow direction, since: where σ is a scalar. The (u × B) = |u||B| · sin(θ ) term is maximised relative to u and B magnitudes when θ = 90 o . There is no singularly clear B-field orientation, however, to achieve this. Since the dipole configuration of B has field lines which vary gradually in orientation, of more significant consideration is that the complex flow field exhibits velocities in many directions at different electrical conductivities. Moreover, the interaction between the flow-field and the MHD forcing affects the cumulative orientation of conductive flow velocity with the magnetic field lines, and so the base flow (without magnetic field) cannot entirely guide predictions of the resultant magnetically influenced flow. These coupled effects, among others, mean that the magnetically affected flow field can behave in counter-intuitive and unpredictable waysas demonstrated through the cases simulated in this study.

VI. CONTROL SURFACE REPLACEMENT -AN OUTLOOK
In light of the complexity of the flow field coupled with MHD forcing effects, design for precise equivalence between magnetic interaction and mechanical surface actuation is a precarious task. However, an important result which can be drawn from these studies, is the identification of classes of configurations which produce morphology equivalence in terms of the flow structure classification.
Since the conditions leading the CLASS A flow structure have been extensively examined, here we demonstrate the two identified classes of cases where CLASS C and CLASS B flow structures emerge: drawing equivalence between surface inclination and magnetic field activation.
Firstly, we introduce the additional case where the second surface angle is declined. The base flow case (no magnetic field) is depicted in FIG. 25, and demonstrates the formation of a CLASS C flow structure: a sonic line (M=1) contour is overlaid, demonstrating that only the boundary layer and separation region are subsonic. The remaining detached shock layer is fully supersonic, downstream of the shock inflection point. The numerical Schlieren images shown in FIG. 26 therefore show two classes of equivalence. The declined second surface angle of θ 2 = 50 o and the magnetic dipole configuration of Z dist = 10 mm produce flow topologies which can both be identified as CLASS C. The magnetic field configuration could be considered to exhibit an MHD diminishing effect, with the same topology of a declined actuation surface which similarly produces a thin detached shock layer.
The second Schlieren compares the computed flow field for the B max = 2.25 T, Z dist = 30 mm case from the dipole location parameter study, with the base flow (no magnetic field) of the θ 2 = 65 o case. As could be expected, both actions (magnetic and geometric) serve to augment the detached shock position further upstream. The particularly interesting result, is that both conditions produce a morphological adaptation, resulting in the CLASS B flow structure. In both cases, the flow deceleration induces the formation of a second triple point (with second jet stream exhibiting a cascade of internal reflected shocks -seen weakly in the numerical Schlieren).
Whilst the same class of flow topology can be identified between cases, finer internal differences remain. The size of the separation regions and extent of the detached shock above the vehicle corner, differ between cases. Perfect matching would not be expected given the notable differences in the composition of forcing dynamics. However, identification of equivalent classes of structures, is very important, since it provides a basis for quantitative analysis (regarding enhancement effects), which otherwise cannot be translated directly between flow structure classes.

VII. CONCLUSION
Numerical simulation capabilities have been developed for MHD affected, hypersonic flows over 2D axisymmetric nonsimple geometries. Producing quantitatively accurate results for this flow type is a challenge due to the sensitivity of feature formation and real gas modelling of the weakly ionised plasma, including the accurate computation of electrical transport properties. The accurate results (as matched with experiment) presented in this work can be attributed to the combination of: effective high resolution numerical methods for strong shock-boundary interactions over arbitrarily complex embedded geometries, and the direct and accurate prediction of electrical conductivity via the extended 19 species air plasma EoS (plasma19X).
Application of the numerical model to studies of the hypersonic double cone configuration, with applied magnetic field, reveal detailed -and at times counter-intuitive -flow field effects.
A classification system is introduced to aid the analysis of emergent flow topologies. The numerical study of magnetic dipole centre location finds that the flow field is sensitive to magnetic field configuration and geometric variations. The study identifies two classes of effects: (1) differences in magnitudes of MHD enhancement effect, and (2) alteration of the flow topology. One such counter-intuitive result is the reduced MHD enhancement effect for increased second conical surface angle. The detailed flow physics for each topological flow structure is described, revealing complex and coupled FIG. 26. Numerical Schlieren (density gradient) depicting two classes of emergent flow structures: top) CLASS C ,and bottom) CLASS B. In each case an equivalence is drawn between: left) magnetic actuation effect, with contours of |B| depicted in grey, and right) mechanical surface actuation via varied θ 2 angles of the geometry. mechanisms of action.
Critically, this paper demonstrates how, for hypersonic flows with complex shock interactions, the MHD affected flow is not only augmented in terms of shock structure position, but may exhibit topological adaptations to the flow structure, as compared with the base case.
In terms of the surface actuation concept, classes of morphological equivalence are identified between the magnetic interaction effect and mechanical control surface actuation. This is of technological significance for emerging magnetic actuation technologies.