Hyperbolic Conduction: A Fast, Physical Conduction Model Implemented in Smoothed Particle Hydrodynamics

We present the first implementation of hyperbolic thermal conduction in smoothed particle hydrodynamics (SPH). Hyperbolic conduction is a physically-motivated alternative to traditional, parabolic conduction. It incorporates a relaxation time, which ensures that heat propagates no faster than a physical signal speed. This allows for larger, Courant like, time steps for explicit schemes. Numerical solutions of the hyperbolic conduction equations require added dissipation to remain stable at discontinuities and we present a novel scheme for this. Test cases include a simple step, the Sod shock tube, the Sedov-Taylor blast, and a super bubble. We demonstrate how longer relaxation times limit conduction, recovering the purely hydrodynamical results, while short relaxation times converge on the parabolic conduction result. We demonstrate that our scheme is stable with explicit Courant-like time steps and can be orders of magnitude faster than explicit parabolic conduction, depending on the application.


INTRODUCTION
The physics of thermal conduction and diffusion are important across a huge range of scales in astrophysics.This includes modelling hot gas in galaxies and clusters (Chandran & Cowley 1998), superbubbles (Weaver et al. 1977), stellar coronae (Gudiksen & Nordlund 2002) and cosmic rays (Snodin et al. 2006).The traditional approach is parabolic conduction and diffusion equations, also called Fickian diffusion, named for Fick's law of diffusion (Fick 1855).
Parabolic conduction allows information to travel instantly, which is un-physical (Axford 1965).A common ad hoc solution is to limit the flux to be less than an estimated maximum possible or saturated flux (Cowie & McKee 1977).The saturated flux is typically set to a fraction of the heat carried if all the electrons moved in the same direction at a typical thermal electron speed relative to the gas.
With parabolic conduction, it is assumed that the heat flux is a function of local gradients in temperature and density only.In practice, the flux takes a finite relaxation time to develop.A more physically appealing approach is to begin with the Boltzmann equations and include this relaxation time (Axford 1965;Gombosi et al. 1993).Gombosi et al. show that the parabolic heat equation is a first order approximation to the Boltzmann equation.Their second order approximation is a set of hyperbolic equations that effectively include a finite propagation speed.Therefore, if we switch to the hyperbolic conduction equations (sometimes referred to as Telegraph equations or Non-Fickian diffusion), we have a model that respects information propagation.This is equivalent to having a built-in, physical limit on the flux and thus avoids the need to apply saturation limits.
The appropriate energy propagation speed depends on the system of interest.For neutral gas, the speed is close to the sound speed.In the case of a plasma, electrons move faster than ions and carry most of the flux, and thus the speed should be close to that of the electrons (where the typical electron speed is roughly an order of magnitude larger than the sound speed,   ) (Spitzer 1956).Magnetic fields make the propagation directionally dependent (Jubelgas et al. 2004).However, where those same magnetic fields are tangled on small scales, the net effect is to dramatically lower the effective transport speed.In general, we expect a transport speed related to the sound speed by a factor dependent on the physical scenario.When the associated relaxation time is short compared to other physics, the precise factor will not strongly affect the behaviour.In the extreme case of very short relaxation times, the behaviour becomes essentially equivalent to parabolic conduction.
An key advantage of hyperbolic conduction is that it is significantly less computationally expensive to use in simulations.Explicit parabolic conduction has very stringent stability limits for its time steps, making it more expensive than other physics (Hanasz & Lesch 2003).Some fluid dynamic codes circumvent the time step issue using implicit methods at the cost of an iterative solution (e.g.Meyer et al. 2012 andDubois &Commerçon 2016).By using hyperbolic conduction, one can maintain a level of mathematical simplicity while allowing time steps similar to those of the hydro solver.Some authors describe hyperbolic conduction as a faster approximation to the parabolic form for numerical work (Rempel 2017).We would argue that it is both faster and more physical.
Prior numerical work has employed hyperbolic conduction in grid codes.Snodin et al. (2006) used it to model cosmic rays in the pencil code.Rempel (2017) and Navarro et al. (2022) have used hyperbolic conduction to model stellar coronae.
In this paper, we will demonstrate a new implementation of hyperbolic conduction, well-suited to smoothed particle hydrodynamics, and incorporated into the gasoline2 code (Wadsley et al. 2017).Hyperbolic equations typically require numerical dissipation for stability near jumps and for this we introduce a modified version of the central scheme in Kurganov & Tadmor (2000) that takes into account the inherent diffusion of the equations.We demonstrate that a Courant-type time step is sufficient for explicit, stable SPH integration.In section three, we present results for standard test problems, including the Sod shock tube and Sedov-Taylor blast.
Given that super bubbles are one of the processes that conduction is used to explore, we also provide this as an example.In the final section, we discuss potential future work.

Hyperbolic Conduction
Parabolic conduction, is governed by the following equations, where  is thermal energy per unit mass,  is the thermal flux,  is mass density, and  is the thermal conductivity.For hyperbolic conduction, we add a time dependant component to equation 1.This leads us to the following, where  is the relaxation time of the system.This is traditionally known as the Telegraph equation and has been verified as a second order approximation to the Boltzmann diffusion equation (Gombosi et al. 1993).Advective terms have been omitted here for simplicity but are included in the full implementation.Because we use a Lagrangian method, advection is included implicitly (Monaghan 1992;Wadsley et al. 2004).
In the cases where  and  are kept constant, these can be rearranged into the traditional heat equation and the hyperbolic version.The parabolic and hyperbolic versions are as follows, (5) Because it is hyperbolic, equation 5 has a characteristic speed that limits propagation.We derive this in appendix B as being To demonstrate how these equations differ, we can look at the general solution for a sample problem. Figure 1 shows the exact solutions (as described in appendix A) for parabolic conduction (dashed lines) and hyperbolic conduction (solid lines) at different times.The initial conditions are a step in temperature from 0.125 to 1 at  = 0. Already, we can see that the hyperbolic transfer hits a virtual wall, limited by  hyp , This virtual wall is more physical than the parabolic process.However, if we do wish to reproduce the parabolic case, we only need to lower our value of  and the two systems become identical as the signal speed increases.
To use these equations properly, we must determine  and .In this first demonstration,  and  are assumed to be constant.However, for some tests that we will perform, we are interested in the behaviour of plasma.For plasma, we use Spitzer conduction to set the , as described in Spitzer (1956).When we solve parabolic conduction explicitly, we use saturation to limit the maximum thermal conduction rates (Cowie & McKee 1977).Because it already limits heat transfer, there is no need to include saturation for the hyperbolic version.Our approach to saturation and Spitzer conductivity are described in Keller et al. (2014).
Once we have ,  can be set using equation 6.The fastest signal speed should be equivalent to the speed of the signal carrying particles (electrons) (Spitzer 1956).We would expect this to be of on the order of the speed of sound.However, Cowie & McKee (1977) note that their saturation equation, and therefore signal propagation, is changed depending on whether an electron's movement aligns with local magnetic fields.Dubois & Commerçon (2016) use a tensor conductivity to allow for anisotropic conduction along field lines.This is beyond the scope of this work but should be feasible.In our large scale simulations, it is commonly assumed that magnetic fields would be tangled below the limit of our resolution.In this case, the isotropic limit is a good representation and that the signal speed reverts to being on the order of the sound speed.With this in mind, we can set  = / 2  , as we would for simulations with weak or no magnetic field.

Stable Time-Steps
Explicit diffusion requires a time-step limit of order 0.25 Δ 2 /, depending on the exact numerical representation of the diffusion term (Hanasz & Lesch 2003).The Δ 2 term is particularly troublesome.It requires that our time step goes down as the square of the resolution (e.g.particle spacing, ℎ).
With hyperbolic conduction, stable timesteps become linked to signal crossing times.This is effectively a standard Courant-Friedrich-Lewy time step criterion (Rempel 2017;Snodin et al. 2006).For hyperbolic conduction, the time-step limits is thus Δ ≤  Δ/ hyp , where  is a factor of order unity ( ∼ 0.5).Keeping in mind our previous derivation, we can also express this as  Δ √︁ /.In cases where  hyp ≈   , Δ =  Δ/  , which is satisfied via the Courant condition that is already used in gasoline2 (Wadsley et al. 2017).
This means that our time step decreases inversely with the particle spacing, thus lengthening our time steps compared to parabolic conduction.The main drawback of using this method is that we now have to keep track of more variables than before and calculate more rates of change at each individual step.However, the larger time steps make up for this.

SPH Implementation
In papers such as Monaghan (1992), parabolic conduction is adapted to SPH as follows, where the right hand side is an SPH approximation for the second derivative in equation 4 and the left hand side is the comoving rate of change of the energy per unit mass,   for particle .In this equation,   ,   and   are the mass, density and conductivity of particle , r   is the vector distance between the two particles and    is the kernel function, (see Monaghan 1992).
For hyperbolic conduction, thermal flux can evolve separately from thermal energy and thus requires its own variable.We adapt equations 2 and 3 separately instead of using equation 5. Our SPH equations for hyperbolic conduction are shown in equations 8 and 9).
In these equations, any variable with of the format ℎ   is the value averaged between particles i and j.The first term on the right hand side of equation 9 is the SPH equivalent of equation 2. The second is a numerical dissipation term inspired by the dissipation term proposed by Kurganov & Tadmor (2000), who note that finite difference-based numerical solutions of hyperbolic equations without such terms generally experience instabilities.These take the form of ringing effects. 1 and  2 are factors slightly less than unity (  1 ∼ 0.5,  2 ∼ 0.1) used to optimize the amount of dissipation.ℎ  is the SPH smoothing length of the particle.If this dissipation were implemented in a grid code, ℎ could be replaced by the grid spacing.
The added dissipation term, as outlined in Kurganov & Tadmor (2000), takes the form of a maximum signal propagation speed multiplied by the resolution, Δ (ℎ in SPH), and the second derivative of thermal energy.A similar process was also proposed earlier by Monaghan (1997) with their introduction of an artificial viscosity term into SPH.
A twist on the usual scenario to which Kurganov & Tadmor applied their dissipation, is the fact that we are modeling a hyperbolic approximation to a diffusion equation.Specifically, as  gets small compared to other timescales in the system, it should closely mimic pure diffusion and finite difference-based models of diffusion equations do not need extra dissipation to remain stable.
When √︁ / > /ℎ, the numerical dissipation becomes greater than the thermal diffusion associated with parabolic form.This issue is strongly apparent in the step test of section 3.1.Therefore, we limit the signal velocity to always be less than or equal to /ℎ.This is equivalent to capping the signal speed when  is small.We demonsrate that this is sufficient to retain stability in section 3.1.

Integration Methods
To integrate these equations, we use the leap-frog symplectic integration method described in the original gasoline and gasoline2 papers (Wadsley et al. 2004(Wadsley et al. , 2017)).gasoline2 uses the kick-driftkick (KDK) form of the leapfrog.Velocities are kicked using calculated forces to a point half a time step ahead.Position is drifted to the end of the time-step using that mid-point velocity.Forces are then updated based on the new positions and velocities are kicked forward to the end of the time-step (Wadsley et al. 2004).If the rate of change of velocity depends on velocity (such as is the case for shocks or viscosity), we need a predicted velocity.This also applies for the flux, .Thus we must predict fluxes to update .This done using / from the previous kick.In the current implementation we are storing a , predicted  and a / for each particle in addition to the other variables already present in gasoline2.
An important aspect of gasoline2 is that each particle has an individual time step which is ultimately rounded down to a power of two fraction of the largest allowed or root time step,  0 , so that Δ =  0 /2  . is referred to as the rung of the particle (Wadsley et al. 2004).Thus each full step is the product of pairwise recursive (KDK) sub-steps on higher rungs.

TEST PROBLEMS
In this section, we will explore a number of tests that we can use to demonstrate the behaviour of hyperbolic conduction generally, and our SPH implementation of it.As noted in Wadsley et al. (2017), lower-dimensional SPH tests do not accurately predict SPH behaviour in 3D simulations.Thus we perform all our tests in 3D.

Step Function
A key test for hyperbolic conduction is a step function.Being diffusive, numerical conduction behaves very well on smooth test problems, so we did not include such tests results here.On the other hand, the presence of an initial step can strongly excite ringing in numerical solutions of hyperbolic equations.In each test in this section, the density, , the thermal conductivity, , and the relaxation time, , were held constant.
The initial condition is of a uniform hot material with an energy per unit mass,  = 1 placed adjacent to a cold one with  = 0.125 with the interface at  = 0. Hydrodynamics are not modeled.The particles are arranged in a glass with a unit volume with 16 3 particles replicated along the x-axis so that the particle spacing is ∼ 0.0625.These are thus 3D tests and are periodic with a period of 1 in the transverse directions.
Figure 2 demonstrates how the solution smoothly varies as we change the relaxation time, , with a fixed thermal conductivity,  = 1.The density is kept at  = 1.The time is set to 1 for each test so that we can see how hyperbolic conduction gradually approaches the parabolic limit as we decrease .In cases with a long relaxation time, the step in temperature propagates to a distance of  hyp  = √︁ 1/ (only in this case as  =  = 1).In cases of small  ≲ 0.1 , the parabolic and hyperbolic tests are hardly distinguishable.
Figure 3 shows a version of this test with the same step initial conditions for .Here, we compare different numerical dissipation expressions used to keep the numerical integration stable.We also vary  and keep  = 1 and  = 0.125 for this test.The output time has been selected so that the step in heat always propagates a unit distance from the origin,  = 1/ hyp = 0.0035 √︁ 1/.Here, we demonstrate the improvement associated with using the minimum of the two characteristic velocities for numeric dispersion expressed in equation 10.
The left panel shows results using /Δ as the characteristic speed.This is a characteristic speed for diffusion within a single resolution element.It is too diffusive when we are far from the parabolic limit (large  in this case or  ≪ ).The middle panel shows a Kurganov & Tadmor (2000) inspired scheme, using the signal speed,  hyp = √︁ /, as the characteristic speed.It is too diffusive as we approach the parabolic limit (small  here or  ≫ ).The right panel shows results using the minimum of these two as the characteristic speed.Taking the minimum speed consistently recovers the exact solution.Although this figure only shows the solution with one value of , tests with other values of  produce equivalent results.
As well as the two depicted tests, we also ran tests with various other values of  and  in order to calibrate the choice of the constants  1 and  2 in equations 10.Simulations with larger values of  (or smaller ) remain stable with larger time steps, and thus a shorter run-time.The run-time is shorter than parabolic conduction, even in many of of the cases where hyperbolic conduction produces the same results as parabolic conduction.This is the exact behaviour that we wished to achieve by introducing hyperbolic conduction.
There should be some degree of caution taken in selecting very small  values when we wish to replicate parabolic conduction.We can create a system that takes longer to run than a parabolic conduction model, but this requires that  is smaller than Δ 2 / (the diffusion time for a single element).From our step tests, as well as from tests in later sections, we find it is not necessary to use  values this small for hyperbolic conduction to be a close approximation of parabolic conduction.

Sod Shock-Tube
The Sod Shock-Tube is a hydrodynamic test that has an initial jump in pressure and density but not velocity.We use initial conditions for the shock tube described in Zamora et al. (2023), although the shock tube is a common test that can be found in other papers including Pakmor et al. (2011), Gottlieb & Groth (1988), and Wadsley et al. (2017).Our shock tube uses values of (1, 0, 1) for the density, velocity and pressure respectively for the left side of the shock.For the right side of the shock, we have initial conditions of (0.125, 0, 0.1).Our particles use the same glass layout described in section 3.1.The  value for both sides of the the shock is 7/5, and the thermal conductivity  is set to 1.
Our results are shown in figure 4, with exact solutions to the Riemann problem shown in black for comparison.The test with  = 0.1 is already a close match to parabolic conduction, while the cases with  ≥ 100 give similar results to a non-conductive (purely hydrodynamical) shock-tube.Because we would like to limit conduction to the sound speed, we also tested the shock-tube in the case where  = / 2  .This is shown in figure 5.In this shock-tube, the sound speed varies between roughly 0.9 and 1.1.
We see promising results right away.As we can see in figure 4, the case where  = 0.1 mimics parabolic conduction.Intuitively, we would expect this because  is much less than .
For the sound-speed derived relaxation time shown in figure 5, reducing the the coefficient on the relaxation time by a factor of 10 causes it to almost perfectly reproduce parabolic conduction.In general, this may be a valuable rule of thumb to use. ∼ 0.1/ 2  ≲  nearly replicates parabolic conduction.There is some additional noise in the low  case but the system remains stable.This gives a signal speed which is about 3 times the sound speed, with correspondingly smaller numerical time steps.
In terms of computational expense, the minimum time step for parabolic conduction is 128 times smaller than the minimum time step for the  = 0.1/ 2  test (being a power 2 due to our time stepping scheme).Every other hyperbolic has runs with higher time steps than this.This is a huge increase in the speed; it can be the difference between a test taking hours and minutes to run.Therefore, even if we would like to perfectly replicate parabolic conduction, a hyperbolic system with low  remains the better option.
Finally, our replication of the non-conductive system shows the level of control hyperbolic conduction provides over the conductive properties of the system, with only small amounts of diffusion compared to the exact solutions.Thus by varying the relaxation time, ,   we can explore disparate physical regimes (and even effectively fully suppress conduction if we wish).

Sedov-Taylor Blast
The Sedov-Taylor Blast is an explosion test dating back to the 1940's (Taylor 1950;Sedov 1959).The test case offers a very strong shock with an analytic solution.Because it is a strong shock, as mentioned in Wadsley et al. (2017), we must select our initial conditions very carefully.For the ideal Sedov-Taylor Blast, for which we have an analytical solution, we would begin with a point explosion and the background temperature would be 0. Our test uses 64 4.2×10 8 kelvin (∼ 10 55 ergs) particles placed in the centre of the box.This box is a glass of 128 3 particles.The surrounding particle energies correspond to 10 −3 kelvin, a number well below the central temperatures.
Figure 6 shows the resulting blasts at three times for the no conduction (top row), parabolic conduction (middle row), and hyperbolic conduction (bottom row).Because this test is meant to represent an explosion within a real plasma, we use the equation for Spitzer Conduction to set  and sound speed to set  (Spitzer 1956;Cowie & McKee 1977).The results near the shock-wave are similar with and without conduction.Inside the bubble, we see a small increase in density at small radii, corresponding to slightly lower temperatures.The amount of affected mass is very small.
Examining the individual particle timesteps, we find that the minimum time step used for the parabolic case was 16 times smaller than the hyperbolic case.For completeness, we note that the conduction free (purely hydrodynamic test) had a minimum time step slightly larger than the conductive runs.Overall, the computational effort for the parabolic run was much larger than the hyperbolic case.
The large temperature contrast of the Sedov-Taylor blast produces strong conduction with high characteristic velocities entering the numerical dissipation expression.This results in more temperature diffusion in the high temperature, low density central region in the hyperbolic case compared to the parabolic case.This lowers the temperatures slightly (as seen in the figure ).
The blast appears to show little change at the shock due to conduction, which is expected as the temperatures are low enough to limit Spitzer conduction.Weaver et al. (1977) first showed that thermal conduction can play an important role in the formation of wind-driven bubbles around star clusters.Conduction allows additional hot mass to be fed into the bubble interior, changing its temperature (Weaver et al. 1977).This model also applies for supernova driven bubbles.It provides a basis for energy-based, sub-grid models for stellar feedback (Keller et al. 2014).

Superbubble
Our test is shown in figure 7. We continuously inject hot particles into the centre of our system to simulate hot supernova ejecta from a star cluster of 30,000 solar masses.The star cluster particle is placed in the centre of a 4x4x4 kpc box with 128 3 particles.This closely follows the direct injection test demonstrated in Keller et al. (2014), which we aim to reproduce with hyperbolic conduction.From left to right, the panels show the non-conductive case, parabolic conduction, and hyperbolic conduction.
There is a difference in the size of the bubbles depending on whether or not we include conduction.Hyperbolic conduction produces similar behaviour to parabolic conduction.The conductive cases have higher levels of hot mass injected back into the bubble.The overall result is more mass inside the bubble and a larger bubble in the conductive cases, as shown in the right two panels of figure 7.Although the projected density inside the bubble is increased in non-conductive case, analysis shows that the mass density inside the bubble is similar in all three cases.The conducting results match those of Keller et al. (2014).We see differences in run time.For parabolic conduction, the minimum time step is 16 times smaller than the minimum timestep for both the hyperbolic and non-conductive run which have the same minimum time step.The increase in speed from switching to hyperbolic conduction is clear.
Our superbubble tests also demonstrate the effect of conduction on an expanding bubble.Vishniac instabilities are expected to grow on the surface of the bubble, just as we see in figure 7 (Vishniac 1983).These instabilities arise from the cool shell of the superbubble interacting with the surrounding medium, They are more developed in cases where the bubble is conductive.

Future Work
Our hyperbolic conduction implementation uses more numerical dissipation near discontinuities than the parabolic version.For example, there is a small temperature decrease in the centre of the Sedov-Taylor blast.Our discontinuity model (and associated parameters) were calibrated to supply enough dissipation to remain stable on tests such as the shock tube (section 3.2).However, with large temperatures, Spitzer conductivity can lead to large characteristic speeds and correspondingly large numerical dissipation.Our numerical dissipation is modeled after the central scheme of Kurganov & Tadmor (2000), who apply it using higher order methods such as linear reconstruction.If we were to use piecewise linear rather than piecewise constant reconstruction we can probably lower the dissipation substantially, particularly for non-zero gradients, as suggested by Beck et al. (2016).It would be worth exploring this as an addition to gasoline2 in the future.
So far, we have assumed that magnetic fields play no direct role in the conductive process.As noted previously, unresolved, tangled magnetic fields can effectively make the diffusion isotropic.More generally, if we wish to use hyperbolic conduction with magnetic fields, we should use anisotropic conduction similar to that found in Dubois & Commerçon (2016).This also applies when modeling cosmic ray diffusion (Snodin et al. 2006).
A key driver for this work has been to facilitate exploration of the role of superbubbles as a mode of energetic stellar feedback (from winds and supernovae).The hot mass in the bubble and its temperature are critically dependent on conduction.Once we move beyond idealized superbubbles (Weaver et al. 1977) to the real ISM, we need to consider a clumpy medium, magnetic fields affecting conduction, unresolved cooling and turbulent conduction.These effects were explored in one dimension by El-Badry et al. (2019).This is a promising potential application of our research.For example, Keller et al. (2014) developed a sub-grid two-phase model based on idealized superbubbles which could benefit from an update that incorporates these non-ideal effects.

Summary
We have argued that hyperbolic conduction is a viable, and physically more appealing, alternative to using parabolic conduction.We have also demonstrated that it works in practice with a robust SPH implementation.
Hyperbolic conduction has the advantage that we can limit the propagation of information to the relevant physical signal speed.Thus changes cannot propagate instantly to any distance which is the primary unphysical aspect to the parabolic form.A related benefit is that we do not need to contemplate applying a saturation cap to the allowed fluxes.An advantage for simulations is that numerical stability only requires a simple Courant-type time step limit based on the signal speed.This allows us to use larger time steps (similar to those for hydrodynamics) and therefore run our simulations faster.
We implemented hyperbolic conduction within SPH.Traditionally, SPH uses explicit numerical dissipation to handle discontinuities (such as artificial viscosity for jumps in the hydrodynamic equations).We found that a standard prescription based on the hyperbolic signal speed created excess diffusion.We demonstrated that the intrinsic diffusion in the hyperbolic conduction equations allows us to use progressively lower numerical dissipation while remaining numerically stable as the relaxation time becomes shorter and we approach the parabolic regime.
In all cases, the time step limit remains Courant-like, Δ ≲ Δ/ hyp , and our code ran stably without the need for the rather prohibitive Δ ≲ (Δ) 2 / time step associated with parabolic diffusion.We note that hyperbolic conduction brings some computational overheads, such as the need to store and evolve the vector flux quantity.
A hyperbolic conduction signal speed comparable to the sound speed generally gives results that are qualitatively and quantitatively similar to the parabolic case.There are good physical arguments to support the heat transport speed being similar to the sound speed.In many test problems, there are other physical timescales that are much slower than these speeds and the precise choice of speed has limited effects.In these cases, the sound speed is an attractive choice.However, the precise signal speed may be higher (e.g. a plasma assuming no magnetic fields).The user might also desire a very close match to the parabolic result, which generally requires a signal speed 3 times the sound speed or faster.In such cases, the hyperbolic approach is still typically more computationally efficient than an explicit parabolic approach.

Figure 1 .
Figure 1.Analytical solutions for parabolic and hyperbolic conduction at different times.This test begins with a step in temperature from 0.125 to 1 at x=0 with a mass density of 1 and  = 1.Hyperbolic is shown as solid lines whereas parabolic is shown as dashed lines.The top image shows the case where  is large enough to stop heat from propagating past a point, limited by the signal speed,  hyp .

Figure 2 .
Figure 2. Hyperbolic conduction implemented in SPH.Initial conditions are a step in temperature from 0.125 to 1 at  = 0.  = 1,  = 1,  = 1, and  varies as indicated.Solid lines are the exact solution for hyperbolic conduction.The dashed black line is the parabolic solution (this largely overlaps with the solution for  = 0.1).Diamonds are the average numerical solution solved using SPH, binned at the particle spacing.Values from individual particles are plotted in grey.

Figure 3 .
Figure 3. Step in temperature with  = 0.0001,  = 0.125, and varying .The left hand panel uses 0.1/Δ as the characteristic speed.The middle uses 0.5 √︁ /.The right panel takes the minimum of both.The exact solutions are represented by the solid lines.Initial conditions and time are the same as in figure 2. In each case, the time is set to the crossing time for a distance of 1 (e.g. = √︁ / = 0.0035 √︁ 1/).

Figure 4 .
Figure 4.The conductive Sod Shock-Tube using different values of  as well as parabolic conduction at  = 1.Initial values of density, velocity, and pressure are (1, 0, 1) respectively on the left of  = 0 and (0.125, 0, 0.1) on the right.The exact solution to the Riemann problem is shown in black.The top row includes cases that more closely resemble purely parabolic conduction while the bottom row shows cases closer to the non-conductive limit.Values from individual particles are plotted in grey while colored particles are bins of with a physical width on the order of one particle spacing.

Figure 5 .
Figure 5.The Sod Shock-Tube with  set using  = 0.1/ 2  (red),  = / 2  , as well as parabolic conduction.Initial conditions and time are the same as described in 4.Not shown is an additional test with  = 0.01/ 2  which almost perfectly replicates parabolic conduction.

Figure 6 .
Figure6.The Sedov-Taylor blast at 1, 2, and 3 Myr.This test includes a plasma of 4.2 × 10 8 kelvin gas in the centre of a cold box.The top row shows no conduction while the bottom two show parabolic and hyperbolic conduction respectively.Solid lines show the exact answer for a point explosion as described in Taylor (1950);Sedov (1959).Black dashed lines show the initial distribution of temperature for our particles.

Figure 7 .
Figure 7.A super bubble test 50 Myr after the creation of the star cluster.Hot mass is continuously injected from stellar winds in a central 30,000  ⊙ star cluster.The top left test is run with no conduction.The centre test uses parabolic conduction.The rightmost test uses hyperbolic conduction with  = / 2  .The conductive bubbles grow much larger than their non-conductive counterpart and create larger instabilities.