2014 Highlights

Figure 1. A volume rendering of the magnetic field distribution in a simulation of a region of quit Sun.

In the outer layers of the solar interior, vigorous convective motions interact with complex magnetic field structures. Regions of quiet Sun (which, by definition, contain no large-scale magnetic features such as sunspots) are characterized by a spatially intermittent distribution of small-scale magnetic fields. Explaining the origin and evolution of quiet Sun magnetic fields is one of the most important challenges in modern solar physics.

It is believed that a significant fraction of quiet Sun magnetism is generated in the surface layers by a convectively-driven dynamo. To simulate this dynamo process, it is necessary to solve the governing equations of Magneto HydroDynamics (MHD) for a compressible, electrically-conducting gas in three spatial dimensions and time. This is a computationally challenging problem that can only be tackled by using massively parallel simulations on the Dirac High Performance Computing facilities.

Figure 2. The temperature distribution (left) and vertical magnetic field distribution (right) from a simulation of a region of quit Sun.

Using simulations of this type, we have confirmed that a dynamo mechanism could indeed account for the observations. In the absence of magnetic fields, the simulated convective flows are solar-like in many respects, with a “granular” pattern of broad warm upflows, surrounded by a network of cooler, narrow downflows. Furthermore, as suggested by the observations, these granules are organized on a larger (“mesogranular”) scale. When a seed magnetic field is introduced into this system, the convective motions drive a dynamo with a complex, intermittent field distribution that is comparable to that observed in the quiet Sun, with the most persistent magnetic structures accumulating preferentially at the mesogranular boundaries.

Plasma turbulence is believed to be responsible for the fact that the solar wind – the supersonic expansion of the Sun’s corona – retains a high temperature throughout the heliosphere, even when it is expected to cool due to expansion. The solar wind may be heated by releasing magnetic energy to accelerate particles at thin magnetic discontinuities, or current sheets. In some models these current sheets are formed at the same time as the solar wind leaves the Sun, or they may form spontaneously as part of the turbulence seen in the solar wind.

Figure 1. Isosurfaces of current sheet density showing rippling due to ion drift-kink instability (right), and (left) virtual spacecraft data for magnetic field for the crossing of the current sheet shown in the top panel. (From Gingell et al., 2014, arXiv: 1411.4422.)

Using DiRAC we have carried out self-consistent, particle kinetic plasma simulations of the evolution of thin current sheets in a three-dimensional geometry. We find that there are two important instabilities which control energy release. The ion tearing instability releases energy via magnetic reconnection, and the drift-kink instability which, as it name implies, produces a rippled current. The drift-kink instability can reduce the rate of energy release, and lead to thicker current sheets which are less liable to reconnect due to the tearing instability. The competition between these two instabilities can only be studied with three-dimensional simulations, and these require large memory and processor numbers which would not be possible without DiRAC resources.

Our results indicate that the importance of thin current sheets to solar wind heating might be less than predicted by two-dimensional simulations which only capture the ion tearing instability. Energy release at current sheets could be important if they are actively produced at a high enough rate by large scale turbulence. We have also created virtual spacecraft data which can be used to test our simulations against actual data from spacecraft.

Figure 1. Numerical simulation of a gas cloud accreting on to a SMBH binary, performed on the DiRAC-2 Complexity cluster. The simulation followed the evolution of a turbulent, self-gravitating cloud of gas as it falls towards the black holes and forms an accretion disk around the binary.

Almost all galaxies have supermassive black holes (SMBHs) in their centres. These cosmic giants are millions to billions of times more massive than the Sun, and play a key role in shaping the formation and evolution of galaxies. Galaxies grow through collisions and mergers of smaller galaxies, and when pairs of galaxies merge their SMBHs slowly sink towards the centre of the merged galaxy. As they get closer together the “friction” dragging them inwards (actually gravitational encounters with stars) becomes progressively less and less efficient, and this process ultimately stalls when the SMBHs reach separations of around a parsec. However, no pairs of binary SMBHs are seen at these separations (despite many years of searching), suggesting that the process(es) that drive SMBH mergers are much more efficient than theory suggests. This long-standing discrepancy between theory and observations is known as the “last parsec problem”.

The most popular solution to this problem is to appeal to a gaseous accretion disc around the SMBH binary. These discs are commonly observed around single SMBHs (where they appear as quasars), and they allow SMBHs to grow by funnelling material down towards the event horizon. However, it is well-known that if the binary and disc are aligned (i.e., they orbit in the same plane and rotate in the same direction), then the disc cannot shrink the binary’s orbit quickly enough. The Theoretical Astrophysics Group in Leicester has recently pioneered the idea that misaligned accretion discs may be the key mechanism driving SMBH mergers. The interaction between an accretion disc and a misaligned SMBH binary is very complicated, and can only be studied using large computer simulations. The power of DiRAC-2 allowed us to build detailed simulations of how accretion discs form and evolve around SMBH binaries (see Figure 1, adapted from Dunhill et al 2014. Animations from these simulations are available here). These are the first models to study the disc formation process self-consistently (previous studies have simply assumed that a disc is present), and to look in detail at how the gravity of the disc shapes the evolution of the binary. As with previous studies, we find that prograde discs (i.e., disc which rotate in the same direction as the binary) do not shrink the binary’s orbit efficiently. However, when discs form in retrograde configurations (i.e., with the gas orbiting in the opposite direction to the binary) their interaction with the binary is much more efficient. Our simulations demonstrated that repeated retrograde accretion events can drive mergers between SMBH binaries on short time-scales, and this process offers an elegant potential solution to the last parsec problem.

Novel strong interactions might provide an explanation for the mechanism of electroweak symmetry breaking. In order for a strong dynamics to be able to explain the electroweak breaking, the theory must be near the onset of the conformal window and possess an anomalous dimension of the chiral condensate of order one. Furthermore, after the experimental observation of a light Higgs and no unexpected particle, a crucial requirement for a Beyond the Standard Model (BSM) theory is the existence of a light scalar in the spectrum.

The a priori interesting theories can be divided into two classes: infrared conformal theories, which are inside the conformal window and are scale-invariant in the limit of infinite distance, and walking theories, showing confinement at large distance, but behaving as infrared conformal theories for a wide range of intermediate scales. Theories in both categories can be used to construct BSM models without violating the existing constraints from precision measurements.

Figure 1. Selected spectral states of the theory in lattice units.

Our investigations of the SU(2) with one adjoint fermion model show that this theory displays crucial dynamical features attributed to phenomenologically viable models of electroweak symmetry breaking (see figure 1 of the spectrum from arXiv:1311.4155).

Indeed simulations performed for values of the couplings in the region connected with the continuum indicate that the system is close to the edge of the conformal window, with an anomalous dimension close to one and a scalar particle much lighter than the other particles in the spectrum.

At the present stage this is the only system showing such intriguing features. Moreover, our calculation provides the first evidence that concrete models with the non-perturbative features mandated by experimental constraints do exist.

Figure 1. The interaction between a muon and a photon is complicated by the production of virtual particles.

A vacuum is never completely empty space but teems with particles that are created fleetingly by quantum fluctuations in energy and then disappear. The heavier the particle, the less time its brief existence can last. In that time, however, the particle can interact with the much lighter particles in our everyday world and thereby reveal its existence through tiny discrepancies in the properties of these particles from that expected in the Standard Model. The magnetic moment of the muon shows such a discrepancy, a tantalizing 25±9 parts in 10¹⁰. To pin this down more accurately a new experiment, E989 at Fermilab near Chicago, will reduce the experimental uncertainty by a factor of 4. An improved theoretical uncertainty from the Standard Model to match this needs lattice QCD calculations. This year HPQCD has developed a new method using DiRAC that will enable the accuracy needed to be achieved.

Figure 1. Our results for the strange quark HVP contribution to aμ as a function of lattice spacing.

The muon is a heavier cousin of the electron with the same electric charge and spin, so that it carries a magnetic moment. The ‘gyromagnetic ratio’, g, measures the ratio of magnetic moment to spin (in units of e/2m). Naively, g=2 but this ignores the interactions of the muon with the cloud of virtual particles discussed above. The anomalous magnetic moment, a_μ = (g – 2)/2, can be measured very accurately from the difference between the circulation and spin precession frequencies as muons fly round a ring transverse to a magnetic field.

The biggest effect in the Standard Model result for a_μ comes from Quantum Electrodynamics where very high order calculations can be done. The largest theoretical uncertainty is from Quantum Chromodynamics and the process of Fig. 1 in which a muon radiates a photon that briefly becomes a quark-antiquark pair. This ‘hadronic vacuum polarisation’ (HVP) contribution needs to be calculated to better than 0.5% in lattice QCD to reduce the theoretical uncertainty. Our new method (arXiv:1403.1778, Phys. Rev D89:114501) makes this possible, having achieved 1% for the first time for the strange quark HVP contribution. Our method has low systematic errors and DiRAC enabled high statistics on realistic gluon field configurations. The up/down quark calculation is now underway. Errors will be reduced by a further factor of 3 as we increase statistics by a factor of 10 in 2015.

Light baryons such as the familiar proton and neutron make up most of the matter in the visible universe. These baryons have constituents of three types of ‘flavour’ of quarks: up, down and strange or u, d and s. Particles group together in families or multiplets, the octet multiplet is shown on the left-hand side of the figure.

The mass splitting between the proton and neutron is a very delicately balanced quantity, partly caused by the mass difference between the u and d quarks, and partly by QED effects. Small changes in this mass difference would have profound effects on the way the universe looks today. The initial balance between hydrogen and helium, established in the first half-hour after the Big Bang, depends strongly on the neutron half-life, and so on the p-n mass splitting. Later, the production of carbon and oxygen in stars also depends strongly on the proton-neutron splitting.

Figure 1. The octet multiplet (left) and a compilation of recent lattice determinations of baryon mass splittings (right).

The strong force binding the quarks together is described by quantum chromo-dynamics or QCD, which has to be numerically simulated (via lattice QCD). A compilation of recent lattice determinations of baryon mass splittings is given in the right panel of the figure. In particular there is also mass splitting between the particles at the centre of the octet the Sigma and Lambda. This is more complicated case as the states mix. We have now extended our previous results to include this case and determined the mixing angle and mass splitting. While the effects of mixing on the masses are very small (second order in the angle), it can be much larger for decay rates.

While the main force between the quarks and antiquarks comes from QCD there is also a contribution from the electromagnetic force, QED, which is usually left out of lattice calculations. We are also doing calculations with both forces included, to see how important the effects of QED are.

Firm observational evidence indicates that supermassive black holes are present at the centre of the majority of galaxies already from very early cosmic times all the way to the present day Universe. Feedback from these black holes in the form of large scale outflows is believed to be one of the key ingredients shaping the evolution and properties of galaxies in our Universe.

Recent observations (Cicone et al. 2014) have detected the emission from spatially extended cold gas around a bright quasar that existed when the Universe was less than 10% of its current age. This cold gas is moving at a very high velocity of the order of 1000km/s and has been detected up to 90,000 light years away from the rapidly accreting black hole that is powering the observed quasar. While this high velocity gas has been interpreted as the signature of a powerful black hole outflow, this observation is in tension with simple theoretical expectations which suggest that while rapidly moving, gas should be very hot instead of cold.

Figure 1. A galaxy hosting a black hole. The hot and rapidly moving outflow originating from the black hole is shown with a black contour, cold gas pockets containing up to a billion Solar masses and moving together with the hot outflow are depicted with grey pixels, and inflowing cosmic web filaments are illustrated in orange hues.

Using the moving-mesh code AREPO, researches at IoA/KICC, Cambridge have performed some of the most realistic cosmological simulations that follow black hole growth and feedback in the early Universe. In simulations, a galaxy hosting a bright quasar undergoes extremely high levels of star formation as it is located at the intersection of cosmic web filaments which deliver copious amount of gas to it. This eventually drives a “starburst phase” where a large number of stars explode as supernovae and their joint effect triggers a galaxy-scale wind. Thus, as the even more powerful black hole-driven outflow leaves the central region of galaxy and propagates outwards, it sweeps over the gas that has been polluted by metals, initially produced within stars and subsequently released by supernova blast waves. These “pockets” of metal enriched gas cause part of the hot quasar-driven outflow to cool to low temperatures. This is illustrated in the figure above where the hot and rapidly moving outflow originating from the black hole is shown with a black contour, cold gas pockets containing up to a billion Solar masses and moving together with the hot outflow are depicted with grey pixels, and inflowing cosmic web filaments are illustrated in orange hues. The black hole and its host galaxy are situated at the centre of the image.

Thanks to the simulations performed at the DiRAC-2 facility in Cambridge, researchers at IoA/KICC have found that self-consistent modelling of both supernovae and black hole outflows lead to the formation of cold gas pockets moving with 1000 km/s which are spatially extended over 90,000 light years (30kpc) in remarkable agreement with the puzzling observations of Cicone et al. 2014.

Figure 1. The image above is a slice through the simulation volume, with the intergalactic gas colour coded from blue to green to red with increasing temperature. The inset images show details in a small region of the simulation centered on a spiral galaxy seen face-on.

The EAGLE simulation project is a Virgo flagship program aimed at creating hi-fidelity hydrodynamic simulations of the galaxy population. It has been a major investment of intellectual effort, manpower and computing resources, with game-changing results, opening a window on the role and physics of baryons in the Universe.

EAGLE uses the latest developments in SPH techniques, with state of the art modelling of sub-grid processes such as star formation, feedback and black hole accretion. Unlike previous projects, our feedback scheme does not involve “de-coupling” feedback particles from hydrodynamic forces, or artificially suppressing cooling. The largest simulation, with a volume of 10⁶ Mpc contains 1000s of well-resolved galaxies similar to the Milky Way and100s of galaxy groups and a handful of galaxy clusters. This huge statistical sample allows us to understand the variety of galaxy formation histories, the impact of galaxy environment and to study rare objects such as highly star forming sub-mm galaxies.

Obviously, a simulation of this size is of no use unless it is able to accurately reproduce the properties of the local galaxy population. EAGLE has been spectacularly successful: for example, the model correctly predicts the evolution of the galaxy mass function and the rising trend of galaxy specific star formation rates (Schaye et al. 2014 at z=0, Furlong et al. 2014 showing predictions at higher z). Numerous papers are in advanced stages of preparation comparing additional galaxy properties, such as HI mass functions (Bahe et al), galaxy colour distributions and luminosity functions (Trayford et al), OVI absorbers (Rahmati et al), AGN luminosity functions (Rosas-Guevara et al). The fidelity of the simulation also makes it a powerful tool for understanding the biases in observational surveys.

Having shown that the simulation’s fidelity is excellent, we can now exploit it as a tool for understanding the physical mechanisms by which galaxies form and evolve. The image above is a slice through the simulation volume, with the intergalactic gas colour coded from blue to green to red with increasing temperature. The inset images show details in a small region of the simulation centred on a spiral galaxy seen face-on.

There are four fundamental forces that describe all known interactions in the universe: gravity; electromagnetism; the weak interaction; and the strong interaction, which is the topic of our research. The strong force binds quarks into hadrons such as protons and neutrons, which in turn form the nuclei of atoms, making up more than 99% of all the known matter in the universe.

Normally, the strong interaction binds quarks so tightly that they never escape from within hadrons. However, at temperatures greater than the deconfining temperature, T_C (which is several trillion Celsius!) it undergoes a substantial change in character, becoming considerably weaker, and quarks become “free”. This new phase of matter is called the “quark-gluon plasma” (QGP) and is presumed to have occurred in the first few microseconds after the Big Bang. Despite the QGP being crucial to the development of the Universe (it was born in this state!) it is poorly understood.

Physicists re-create a mini-version of the QGP by colliding large nuclei (such as gold or lead) in particle accelerators at virtually the speed of light. These experiments have been performed in the Large Hadron Collider at CERN and at the Brookhaven National Laboratory in the USA. The region of QGP created in these experiments is incredibly small – just the size of the colliding nuclei. The QGP “fireball” formed expands and cools very rapidly, quickly returning to normal matter where quarks are tightly bound inside hadrons.

Most hadrons melt in the QGP meaning that the quarks that were inside them become free, like satellites escaping the earth’s gravity. However it was proposed 30 years ago that some, heavy “mesons” may remain bound up to quite high temperatures. Our collaboration has calculated the spectral functions of these bottomonium mesons (comprised of two bottom quarks) and we show the results above. Each pane has the results from two nearby temperatures, with the coldest pair in the top-left, and the hottest bottom-right. The strong peak at around 9.5 GeV shows a bound state which gets weaker as the temperatures increase, but remains in existence up to around 2T_C. This confirms the decades-old hypothesis that these mesons do not melt as soon as the QGP is formed, and that measuring their features is an effective way of determining the system’s temperature, i.e. they can be used as a “thermometer” of the QGP.

COSMOS consortium researchers have been exploiting the DiRAC HPC Facilities to make progress towards ambitious milestones in five key inter-related areas: (i) extreme universe, (ii) cosmic microwave sky, (iii) dark energy and (iv) galaxy formation and (v) black holes.

Figure 1. CMB temperature (left) and E-mode polarization (right) bispectrum reconstruction obtained from the Planck Full Mission Planck data using the MODAL estimator on COSMOS. Isosurfaces are shown both positive (red) and negative (blue).

Planck Satellite Science with Polarization:

The COSMOS Consortium flagship project to analyse and interpret data from ESA’s Planck Satellite has gone forward, building on the success of the 2013 First Data Release, to analyse the Planck Full Mission data with CMB polarization in results announced in December 2014. DiRAC resources were used to derive accurate estimates for cosmological parameters, showing interesting deviations from earlier WMAP results, checking consistency with temperature using the polarization data. Planck data was used to impose world-leading limits on cosmic non-Gaussianity that substantially constrain inflationary Universe paradigms. The three-point correlator (bispectrum) was evaluated to high precision for the first time in both temperature and polarization (and mixed bispectra) with a wide-ranging model survey undertaken, which yielded tight constraints on standard local, equilateral and orthogonal non-Gaussianity, as well as a measurement of the lensing bispectrum.

Extreme Universe – Blackholes:

Consortium members continue with world-leading black hole simulations, developing new methods to study black-hole collisions – events so violent their output eclipses the entire electromagnetic universe. Our simulations on Cosmos demonstrate that up to 50% of the total energy can be converted into gravitational waves in such collisions, supporting present assumptions in the analysis of collision search experiments at the LHC. We have also been preparing gravitational wave templates for the Advanced LIGO experiment which will begin operating in 2015. In addition, we have developed and tested a numerical relativity code GRChombo, the first full adaptive mesh refinement code which can be used to simulate black hole collisions and gravitational wave generation in early universe cosmology.

Galaxy Surveys, Dark Energy and Dark Matter:

Consortium resources have been deployed to analyse galaxy surveys through N-body simulations to better understand the properties of the universe and its perturbations. Much of the focus has been on preparing for The Dark Energy Survey (DES), as well as ambitious future surveys with DESI, Euclid and SKA. Power spectrum measurements from Baryon Oscillation Spectroscopic Survey (BOSS) data and mock catalogues corresponding to the Data Releases 10-12 samples were performed using COSMOS.

We have pursued the first ever Quantum Chromodynamics (QCD), or strong force, lattice simulations of quarks and gluons with continuum limit results from chiral lattice quarks at their physical masses. Only left handed quarks couple to the electroweak forces. Our chiral approach reproduces the handedness symmetry of continuum QCD and is particularly good for determining the rates of weak decay processes involving left handed quark operators. This preference for left handed over right handed shows that the weak forces surprisingly break reflection (or parity) symmetry. Our results confirm the strong force as the confining theory of hadrons (particles made from quarks). The mass spectrum and pion decay constant agree with experimental determinations at the fraction of a percent level, and allow determination of the V_us quark flavor mixing parameter of the standard model. These and the following achievements are crucial ingredients in interpreting the results from the large experimental facilities (e.g. at CERN).

Figure 1. The figure above left shows the fit to simulation data for the decay constant from arXiv: 1411.7017 and right the action density for one of the gauge fields generated on the DiRAC facility.

Quark Masses:

We have determined the up/down and strange quark masses – fundamental constants of nature – to a fraction of a percent.

Neutral Kaon oscillations:

In the Standard Model (SM) of particle physics, weak interactions allow quarks to change from one type or flavour to another in a highly constrained pattern. Quark flavour-changing processes with no accompanying change in electric charge are very rare and present a small background against which to search for new physics. The SM allows electrically-neutral kaons to oscillate into their own antiparticles in a way which breaks the combined symmetry of particle-antiparticle exchange (charge) and mirror reflection (parity), known as CP symmetry. CP violation is essential to enable the dominance of matter over antimatter in the Universe and was first discovered in the neutral kaon system. We have calculated the neutral Kaon oscillation amplitude BK to 0.3% accuracy.

Algorithms:

Vigorous activity by UKQCD has leveraged STFC investment in DiRAC with around a 200x speedup from code and algorithmic improvements:

• UKQCD scientists worked for five years with IBM research to help design the BlueGene/Q chip (being responsible for 8% of the chip area);
• UKQCD have developed new multi-grid algorithms that adaptively accelerate the calculation of quark propagators;
• A volume averaging technique developed by RBC-UKQCD, known as AMA, has lowered our statistical errors to around the 0.2% level.

This project has received the bulk of our effort over the past 12 months. We have been running 127 million particle SPH runs of a spiral shock to study the formation of molecular clouds and establish the first self-consistent initial conditions for simulations of star formation. These simulations follow the formation of molecular clouds along a 4kpc length of a spiral arm with a total mass of 20 million solar masses and can resolve the formation of self-gravitating structures down to 10 solar masses. We have been studying the importance of the thermal physics in triggering star formation as well as several tests to ensure that the dynamics of the shocked regions are properly captured. A companion study using grid-based simulations has been completed which investigates the dynamics of single clouds entering into a spiral arm. This study showed that the gravitational potential, plus shock induced hydrodynamical instabilities, can drive turbulence into the cloud. This re-enforces the conclusions from our Smoothed Particle Hydrodynamics simulations.

We are now developing full galactic models of the ISM including live potentials. This work will form a significant part of our simulation work in 2015.

The Galactic Centre:

We have concluded studies on star formation in the Galactic centre and have expanded our scope to study the dynamics of the Central Molecular Zone including the “Moilinari” ring. These simulations show how a single cloud can reproduce the overall features. A detailed comparison with observations is ongoing. We have also used larger scale simulations of the inner 500pc of the galaxy to study how the non-axisymmetric structures combine with infalling gas interact, in order to probe the relationship between the CMZ and the inner star formation activity.