Neutron stars, and magnetars in particular, are known to host the strongest magnetic fields in the Universe. The origin of these strong fields is a matter of controversy. In this preliminary work, via numerical simulations, we study, for the first time in non-ideal general relativistic magnetohydrodynamic (GRMHD) regime, the growth of the magnetic field due to the action of the mean-field dynamo due to sub-scale, unresolved turbulence. The dynamo process, combined with the differential rotation of the (proto-)star, is able to produce an exponential growth of any initial magnetic seed field up to the values required to explain the observations. By varying the dynamo coefficient we obtain different growth rates. We find a quasi-linear dependence of the growth rates on the intensity of the dynamo. Furthermore, the time interval in which exponential growth occurs and the growth rates also seems to depend on the initial configuration of the magnetic field.

Among the possible extensions of General Relativity that have been put forward in order to address some long standing issues in our understanding of the Universe, Scalar-Tensor Theories have received a lot of attention for their simplicity. Interestingly, some of these predict a potentially observable non-linear phenomenon, known as \textit{spontaneous scalarisation}, in the presence of highly compact matter distributions, like the case of neutron stars. Neutron stars are ideal laboratories to investigate the properties of matter under extreme conditions, and in particular they are known to harbour the strongest magnetic fields in the Universe. Here, for the first time, we present a detailed study of magnetised neutron stars in Scalar-Tensor Theories. First, we show that the formalism developed for the study of magnetised neutron stars in General Relativity, based on the \textit{eXtended Conformally Flat Condition}, can easily be extended in the presence of a non-minimally coupled scalar field, retaining many of its numerical advantages. We then carry out a study of the parameter space considering the two extreme geometries of purely toroidal and purely poloidal magnetic fields, varying both the strength of the magnetic field and the intensity of scalarisation. We compare our results with magnetised general-relativistic solutions and un-magnetised scalarised solutions, showing how the mutual interplay between magnetic and scalar fields affect the magnetic and the scalarisation properties of neutron stars. In particular, we focus our discussion on magnetic deformability, maximum mass and range of scalarisation.

Pulsars out of their parent SNR directly interact with the ISM producing so called Bow-Shock Pulsar Wind Nebulae, the relativistic equivalents of the heliosphere/heliotail system. These have been directly observed from Radio to X-ray, and are found also associated to TeV halos, with a large variety of morphologies. They offer a unique environment where the pulsar wind can be studied by modelling its interaction with the surrounding ambient medium, in a fashion that is different/complementary from the canonical Plerions. These systems have also been suggested as the possible origin of the positron excess detected by AMS and PAMELA, in contrast to dark matter. I will present results from 3D Relativistic MHD simulations of such nebulae. On top of these simulations we computed the expected emission signatures, the properties of high energy particle escape, the role of current sheets in channeling cosmic rays, the level of turbulence and magnetic amplification, and how they depend on the wind structure and magnetisation.

The constrained transport (CT) method reflects the state of the art numerical technique for preserving the divergence-free condition of magnetic field to machine accuracy in multi-dimensional MHD simulations performed with Godunov-type, or upwind, conservative codes. The evolution of the different magnetic field components, located at zone interfaces using a staggered representation, is achieved by calculating the electric field components at cell edges, in a way that has to be consistent with the Riemann solver used for the update of cell-centered fluid quantities at interfaces. Albeit several approaches have been undertaken, the purpose of this work is, on the one hand, to compare existing methods in terms of robustness and accuracy and, on the other, to extend the upwind contrained transport (UCT) method by Londrillo & Del Zanna (2004) and Del Zanna et al. (2007) for the systematic construction of new averaging schemes using the information available from 1D Riemann solvers. Our results are presented here in the context of second-order schemes for classical MHD, but they can be easily generalized to higher than second order schemes, either based on finite volumes or finite differences, and to other physical systems retaining the same structure of the equations, such as that of relativistic or general relativistic MHD.

The recent imaging of the M87 black hole at millimeter wavelengths by the Event Horizon Telescope (EHT) collaboration has triggered a renewed interest in numerical models for the accretion of magnetized plasma in the regime of general relativistic magnetohydrodynamics (GRMHD). Here non-ideal simulations, including both the resistive effects and, above all, the mean-field dynamo action due to sub-scale, unresolved turbulence, are applied for the first time to such systems in the fully nonlinear regime. Combined with the differential rotation of the disk, the dynamo process is able to produce an exponential growth of any initial seed magnetic field up to the values required to explain the observations, when the instability tends to saturate even in the absence of artificial quenching effects. Before reaching the final saturation stage we observe a secondary regime of exponential growing, where the magnetic field increases more slowly due to accretion, which is modifying the underlying equilibrium. By varying the dynamo coefficient we obtain different growth rates, though the field seems to saturate at approximately the same level, at least for the limited range of parameters explored here, providing substantial values for the MAD parameter for magnetized accretion. For reasonable values of the central mass density and the commonly employed recipes for synchrotron emission by relativistically hot electrons, our model is able to reproduce naturally the observed flux of Sgr A*, the next target for EHT.

Recent developments in compact object astrophysics, especially the discovery of merging neutron stars by LIGO, the imaging of the black hole in M87 by the Event Horizon Telescope (EHT) and high precision astrometry of the Galactic Center at close to the event horizon scale by the GRAVITY experiment motivate the development of numerical source models that solve the equations of general relativistic magnetohydrodynamics (GRMHD). Here we compare GRMHD solutions for the evolution of a magnetized accretion flow where turbulence is promoted by the magnetorotational instability from a set of nine GRMHD codes: Athena++, BHAC, Cosmos++, ECHO, H-AMR, iharm3D, HARM-Noble, IllinoisGRMHD and KORAL. Agreement between the codes improves as resolution increases, as measured by a consistently applied, specially developed set of code performance metrics. We conclude that the community of GRMHD codes is mature, capable, and consistent on these test problems.

We describe a novel Godunov-type numerical method for solving the equations of resistive relativistic magnetohydrodynamics. In the proposed approach, the spatial components of both magnetic and electric fields are located at zone interfaces and are evolved using the constrained transport formalism. Direct application of Stokes' theorem to Faraday's and Ampere's laws ensures that the resulting discretization is divergence-free for the magnetic field and charge-conserving for the electric field. Hydrodynamic variables retain, instead, the usual zone-centred representation commonly adopted in finite-volume schemes. Temporal discretization is based on Runge-Kutta implicit-explicit (IMEX) schemes in order to resolve the temporal scale disparity introduced by the stiff source term in Ampere's law. The implicit step is accomplished by means of an improved and more efficient Newton-Broyden multidimensional root-finding algorithm. The explicit step relies on a multidimensional Riemann solver to compute the line-averaged electric and magnetic fields at zone edges and it employs a one-dimensional Riemann solver at zone interfaces to update zone-centred hydrodynamic quantities. For the latter, we introduce a five-wave solver based on the frozen limit of the relaxation system whereby the solution to the Riemann problem can be decomposed into an outer Maxwell solver and an inner hydrodynamic solver. A number of numerical benchmarks demonstrate that our method is superior in stability and robustness to the more popular charge-conserving divergence cleaning approach where both primary electric and magnetic fields are zone-centered. In addition, the employment of a less diffusive Riemann solver noticeably improves the accuracy of the computations.

This paper briefly reviews a number of fundamental measurements that need to be made in order to characterize turbulence in space plasmas such as the solar wind. It has long been known that many of these quantities require simultaneous multipoint measurements to attain a proper characterization that would reveal the fundamental physics of plasma turbulence. The solar wind is an ideal plasma for such an investigation, and it now appears to be technologically feasible to carry out such an investigation, following the pioneering Cluster and MMS missions. Quantities that need to be measured using multipoint measurements include the two-point, two-time second correlation function of velocity, magnetic field and density, and higher order statistical objects such as third and fourth order structure functions. Some details of these requirements are given here, with a eye towards achieving closure on fundamental questions regarding the cascade rate, spectral anisotropy, characteristic coherent structures, intermittency, and dissipation mechanisms that describe plasma turbuelence, as well as its variability with plasma parameters in the solar wind. The motivation for this discussion is the current planning for a proposed Helioswarm mission that would be designed to make these measurements,leading to breakthrough understanding of the physics of space and astrophysical turbulence.

The exponential amplification of initial seed magnetic fields in relativistic plasmas is a very important topic in astrophysics, from the conditions in the early Universe to the interior of neutron stars. While dynamo action in a turbulent plasma is often invoked, in the last years a novel mechanism of quantum origin has gained increasingly more attention, namely the Chiral Magnetic Effect (CME). This has been recognized in semi-metals and it is most likely at work in the quark-gluon plasma formed in heavy-ion collision experiments, where the highest magnetic fields in nature, up to B~10^18 G, are produced. This effect is expected to survive even at large hydrodynamical/MHD scales and it is based on the chiral anomaly due to an imbalance between left- and right-handed relativistic fermions in the constituent plasma. Such imbalance leads to an electric current parallel to an external magnetic field, which is precisely the same mechanism of an alpha-dynamo action in classical MHD. Here we extend the close parallelism between the chiral and the dynamo effects to relativistic plasmas and we propose a unified, fully covariant formulation of the generalized Ohm's law. Moreover, we derive for the first time the 3+1 general relativistic MHD equations for a chiral plasma both in flat and curved spacetimes, in view of numerical investigation of the CME in compact objects, especially magnetars, or of the interplay among the non-ideal magnetic effects of dynamo, the CME and reconnection.

Magnetic reconnection is believed to be the driver of many explosive phenomena in Astrophysics, from solar to gamma-ray flares in magnetars and in the Crab nebula. However, reconnection rates from classic MHD models are far too slow to explain such observations. Recently, it was realized that when a current sheet gets sufficiently thin, the reconnection rate of the tearing instability becomes "ideal", in the sense that the current sheet destabilizes on the "macroscopic" Alfv\'enic timescales, regardless of the Lundquist number of the plasma. Here we present 2D compressible MHD simulations in the classical, Hall, and relativistic regimes. In particular, the onset of secondary tearing instabilities is investigated within Hall-MHD for the first time. In the frame of relativistic MHD, we summarize the main results from Del Zanna et al. [1]: the relativistic tearing instability is found to be extremely fast, with reconnection rates of the order of the inverse of the light crossing time, as required to explain the high-energy explosive phenomena.

Neutron stars are among the most fascinating astrophysical sources, being characterized by strong gravity, densities about the nuclear one or even above, and huge magnetic fields. Their observational signatures can be extremely diverse across the electromagnetic spectrum, ranging from the periodic and low-frequency signals of radio pulsars, up to the abrupt high-energy gamma-ray flares of magnetars, where energies of ~10^46 erg are released in a few seconds. Fast-rotating and highly magnetized neutron stars are expected to launch powerful relativistic winds, whose interaction with the supernova remnants gives rise to the non-thermal emission of pulsar wind nebulae, which are known cosmic accelerators of electrons and positrons up to PeV energies. In the extreme cases of proto-magnetars (magnetic fields of ~10^15 G and millisecond periods), a similar mechanism is likely to provide a viable engine for the still mysterious gamma-ray bursts. The key ingredient in all these spectacular manifestations of neutron stars is the presence of strong magnetic fields in their constituent plasma. Here we will present recent updates of a couple of state-of-the-art numerical investigations by the high-energy astrophysics group in Arcetri: a comprehensive modeling of the steady-state axisymmetric structure of rotating magnetized neutron stars in general relativity, and dynamical 3-D MHD simulations of relativistic pulsar winds and their associated nebulae.

Motivated by recent observational searches of sub-structure in starless molecular cloud cores, we investigate the evolution of density perturbations on scales smaller than the Jeans length embedded in contracting isothermal clouds, adopting the same formalism developed for the expanding Universe and the solar wind. We find that initially small amplitude, Jeans-stable perturbations (propagating as sound waves in the absence of a magnetic field), are amplified adiabatically during the contraction, approximately conserving the wave action density, until they either become nonlinear and steepen into shocks at a time $t_{\rm nl}$, or become gravitationally unstable when the Jeans length decreases below the scale of the perturbations at a time $t_{\rm gr}$. We evaluate analytically the time $t_{\rm nl}$ at which the perturbations enter the non-linear stage using a Burgers' equation approach, and we verify numerically that this time marks the beginning of the phase of rapid dissipation of the kinetic energy of the perturbations. We then show that for typical values of the rms Mach number in molecular cloud cores, $t_{\rm nl}$ is smaller than $t_{\rm gr}$, and therefore density perturbations likely dissipate before becoming gravitational unstable. Solenoidal modes grow at a faster rate than compressible modes, and may eventually promote fragmentation through the formation of vortical structures.

Geometrically thick tori with constant specific angular momentum have been widely used in the last decades to construct numerical models of accretion flows onto black holes. Such discs are prone to a global non-axisymmetric hydrodynamic instability, known as Papaloizou-Pringle instability (PPI), which can redistribute angular momentum and also lead to an emission of gravitational waves. It is, however, not clear yet how the development of the PPI is affected by the presence of a magnetic field and by the concurrent development of the magnetorotational instability (MRI). We present a numerical analysis using three-dimensional GRMHD simulations of the interplay between the PPI and the MRI considering, for the first time, an analytical magnetized equilibrium solution as initial condition. In the purely hydrodynamic case, the PPI selects as expected the large-scale $m=1$ azimuthal mode as the fastest growing and non-linearly dominant mode. However, when the torus is threaded by a weak toroidal magnetic field, the development of the MRI leads to the suppression of large-scale modes and redistributes power across smaller scales. If the system starts with a significantly excited $m=1$ mode, the PPI can be dominant in a transient phase, before being ultimately quenched by the MRI. Such dynamics may well be important in compact star mergers and tidal disruption events.

The extraordinary energetic activity of magnetars is usually explained in terms of dissipation of a huge internal magnetic field of the order of $10^{15-16}$G. How such a strong magnetic field can originate during the formation of a neutron star is still subject of active research. An important role can be played by fast rotation: if magnetars are born as millisecond rotators dynamo mechanisms may efficiently amplify the magnetic field inherited from the progenitor star during the collapse. In this case, the combination of rapid rotation and strong magnetic field determine the right physical condition not only for the development of a powerful jet driven explosion, manifesting as a gamma ray burst, but also for a copious gravitational waves emission. Strong magnetic fields are indeed able to induce substantial quadrupolar deformations in the star. In this paper we analyze the joint effect of rotation and magnetization on the structure of a polytropic and axisymmetric neutron star, within the ideal magneto-hydrodynamic regime. We will consider either purely toroidal or purely poloidal magnetic field geometries. Through the sampling of a large parameter space, we generalize previous results in literature, inferring new quantitative relations that allow for a parametrization of the induced deformation, that takes into account also the effects due to the stellar compactness and the current distribution. Finally, in the case of purely poloidal field, we also discuss how different prescriptions on the surface charge distribution (a gauge freedom) modify the properties of the surrounding electrosphere and its physical implications.

Pulsar Wind Nebulae (PWNe) constitute an ideal astrophysical environment to test our current understanding of relativistic plasma processes. It is well known that magnetic fields play a crucial role in their dynamics and emission properties. At present, one of the main issues concerns the level of magnetic turbulence present in these systems, which in the absence of space resolved X-ray polarization measures cannot be directly constrained. In this work we investigate, for the first time using simulated synchrotron maps, the effect of a small scale fluctuating component of the magnetic field on the emission properties in X-ray. We illustrate how to include the effects of a turbulent component in standard emission models for PWNe, and which consequences are expected in terms of net emissivity and depolarization, showing that the X-ray surface brightness maps can provide already some rough constraints. We then apply our analysis to the Crab and Vela nebulae and, by comparing our model with Chandra and Vela data, we found that the typical energies in the turbulent component of the magnetic field are about 1.5 to 3 times the one in the ordered field.

Pulsar Wind Nebulae, and the Crab nebula in particular, are the best cosmic laboratories to investigate the dynamics of magnetized relativistic outflows and particle acceleration up to PeV energies. Multidimensional MHD modeling by means of numerical simulations has been very successful at reproducing, to the very finest details, the innermost structure of these synchrotron emitting nebulae, as observed in the X-rays. Therefore, the comparison between the simulated source and observations can be used as a powerful diagnostic tool to probe the physical conditions in pulsar winds, like their composition, magnetization, and degree of anisotropy. However, in spite of the wealth of observations and of the accuracy of current MHD models, the precise mechanisms for magnetic field dissipation and for the acceleration of the non-thermal emitting particles are mysteries still puzzling theorists to date. Here we review the methodologies of the computational approach to the modeling of Pulsar Wind Nebulae, discussing the most relevant results and the recent progresses achieved in this fascinating field of high-energy astrophysics.

In the last decade, the relativistic magnetohydrodynamic (MHD) modelling of pulsar wind nebulae, and of the Crab nebula in particular, has been highly successful, with many of the observed dynamical and emission properties reproduced down to the finest detail. Here, we critically discuss the results of some of the most recent studies: namely the investigation of the origin of the radio emitting particles and the quest for the acceleration sites of particles of different energies along the termination shock, by using wisps motion as a diagnostic tool; the study of the magnetic dissipation process in high magnetization nebulae by means of new long-term three-dimensional simulations of the pulsar wind nebula evolution; the investigation of the relativistic tearing instability in thinning current sheets, leading to fast reconnection events that might be at the origin of the Crab nebula gamma-ray flares.

Magnetic reconnection in laboratory, space and astrophysical plasmas is often invoked to explain explosive energy release and particle acceleration. However, the timescales involved in classical models within the macroscopic MHD regime are far too slow to match the observations. Here we revisit the tearing instability by performing visco-resistive two-dimensional numerical simulations of the evolution of thin current sheets, for a variety of initial configurations and of values of the Lunquist number $S$, up to $10^7$. Results confirm that when the critical aspect ratio of $S^{1/3}$ is reached in the reconnecting current sheets, the instability proceeds on ideal (Alfv\'enic) macroscopic timescales, as required to explain observations. Moreover, the same scaling is seen to apply also to the local, secondary reconnection events triggered during the nonlinear phase of the tearing instability, thus accelerating the cascading process to increasingly smaller spatial and temporal scales. The process appears to be robust, as the predicted scaling is measured both in inviscid simulations and when using a Prandtl number $P=1$ in the viscous regime.

We investigate the possible implications of quark deconfinement on the phenomenology of Long Gamma-Ray Bursts focusing, in particular, on the possibility to describe multiple prompt emission phases in the context of the proto-magnetar model. Starting from numerical models of rotating Hadron Stars and Quark Stars in full general relativity we track the electromagnetic spin-down evolution in both the hadronic and quark phase, linking the two families through conservation of baryon number and angular momentum. We give estimates of the timescales and the energetics involved in the spin-down process deriving, in the relevant spin range, the relation between the initial and the final masses and rotational energies, whenever hadron-quark conversion is possible. We show how the results can be used in relevant astrophysical cases such as the double burst GRB 110709B.

Fast reconnection operating in magnetically dominated plasmas is often invoked in models for magnetar giant flares, for magnetic dissipation in pulsar winds, or to explain the gamma-ray flares observed in the Crab nebula, hence its investigation is of paramount importance in high-energy astrophysics. Here we study, by means of two dimensional numerical simulations, the linear phase and the subsequent nonlinear evolution of the tearing instability within the framework of relativistic resistive magnetohydrodynamics, as appropriate in situations where the Alfven velocity approaches the speed of light. It is found that the linear phase of the instability closely matches the analysis in classical MHD, where the growth rate scales with the Lundquist number S as S^-1/2, with the only exception of an enhanced inertial term due to the thermal and magnetic energy contributions. In addition, when thin current sheets of inverse aspect ratio scaling as S^-1/3 are considered, the so-called "ideal" tearing regime is retrieved, with modes growing independently on S and extremely fast, on only a few light crossing times of the sheet length. The overall growth of fluctuations is seen to solely depend on the value of the background Alfven velocity. In the fully nonlinear stage we observe an inverse cascade towards the fundamental mode, with Petschek-type supersonic jets propagating at the external Alfven speed from the X-point, and a fast reconnection rate at the predicted value R~(ln S)^-1.

Classical MHD reconnection theories, both the stationary Sweet-Parker model and the tearing instability, are known to provide rates which are too slow to explain the observations. However, a recent analysis has shown that there exists a critical threshold on current sheet's thickness, namely a/L~S^(-1/3), beyond which the tearing modes evolve on fast macroscopic Alfvenic timescales, provided the Lunquist number S is high enough, as invariably found in solar and astrophysical plasmas. Therefore, the classical Sweet-Parker scenario, for which the diffusive region scales as a/L~S^(-1/2) and thus can be up to ~100 times thinner than the critical value, is likely to be never realized in nature, as the current sheet itself disrupts in the elongation process. We present here two-dimensional, compressible, resistive MHD simulations, with S ranging from 10^5 to 10^7, that fully confirm the linear analysis. Moreover, we show that a secondary plasmoid instability always occurs when the same critical scaling is reached on the local, smaller scale, leading to a cascading explosive process, reminiscent of the flaring activity.

The phenomenology of the emission of pulsars and magnetars depends dramatically on the structure and properties of their magnetic field. In particular it is believed that the outbursting and flaring activity observed in AXPs and SRGs is strongly related to their internal magnetic field. Recent observations have moreover shown that charges are present in their magnetospheres supporting the idea that their magnetic field is tightly twisted in the vicinity of the star. In principle these objects offer a unique opportunity to investigate physics in a regime beyond what can be obtained in the laboratory. We will discuss the properties of equilibrium models of magnetized neutron stars, and we will show how internal and external currents can be related. These magnetic field configurations will be discussed considering also their stability, relevant for their origin and possibly connected to events like SNe and GRBs. We will also show what kind of deformations they induce in the star, that could lead to emission of gravitational waves. In the case of a twisted magnetosphere we will show how the amount of twist regulates their general topology. A general formalism based on the simultaneous numerical solution of the general relativistic Grad-Shafranov equation and Einstein equations will be presented.

We study the linear and nonlinear evolution of the tearing instability on thin current sheets by means of two-dimensional numerical simulations, within the framework of compressible, resistive magnetohydrodynamics. In particular we analyze the behavior of current sheets whose inverse aspect ratio scales with the Lundquist number $S$ as $S^{-1/3}$. This scaling has been recently recognized to yield the threshold separating fast, ideal reconnection, with an evolution and growth which are independent of $S$ provided this is high enough, as it should be natural having the ideal case as a limit for $S\to\infty$. Our simulations confirm that the tearing instability growth rate can be as fast as $\gamma\approx 0.6\,{\tau_A}^{-1}$, where $\tau_A$ is the ideal Alfv\'enic time set by the macroscopic scales, for our least diffusive case with $S=10^7$. The expected instability dispersion relation and eigenmodes are also retrieved in the linear regime, for the values of $S$ explored here. Moreover, in the nonlinear stage of the simulations we observe secondary events obeying the same critical scaling with $S$, here calculated on the \emph{local}, much smaller lengths, leading to increasingly faster reconnection. These findings strongly support the idea that in a fully dynamic regime, as soon as current sheets develop, thin and reach this critical threshold in their aspect ratio, the tearing mode is able to trigger plasmoid formation and reconnection on the local (ideal) Alfv\'enic timescales, as required to explain the explosive flaring activity often observed in solar and astrophysical plasmas.

We present a new code for performing general-relativistic radiation-hydrodynamics simulations of accretion flows onto black holes. The radiation field is treated in the optically-thick approximation, with the opacity contributed by Thomson scattering and thermal bremsstrahlung. Our analysis is concentrated on a detailed numerical investigation of hot two-dimensional, Bondi-Hoyle accretion flows with various Mach numbers. We find significant differences with respect to purely hydrodynamical evolutions. In particular, once the system relaxes to a radiation-pressure dominated regime, the accretion rates become about two orders of magnitude smaller than in the purely hydrodynamical case, remaining however super-Eddington as are the luminosities. Furthermore, when increasing the Mach number of the inflowing gas, the accretion rates become smaller because of the smaller cross section of the black hole, but the luminosities increase as a result a stronger emission in the shocked regions. Overall, our approach provides the first self-consistent calculation of the Bondi-Hoyle luminosity, most of which is emitted within r~100 M from the black hole, with typical values L/L_Edd ~ 1-7, and corresponding energy efficiencies eta_BH ~ 0.09-0.5. The possibility of computing luminosities self-consistently has also allowed us to compare with the bremsstrahlung luminosity often used in modelling the electromagnetic counterparts to supermassive black-hole binaries, to find that in the optically-thick regime these more crude estimates are about 20 times larger than our radiation-hydrodynamics results.

The Crab Nebula is one of the most efficient accelerators in the Galaxy and the only galactic source showing direct evidence of PeV particles. In spite of this, the physical process behind such effective acceleration is still a deep mystery. While particle acceleration, at least at the highest energies, is commonly thought to occur at the pulsar wind termination shock, the properties of the upstream flow are thought to be non-uniform along the shock surface, and important constraints on the mechanism at work come from exact knowledge of where along this surface particles are being accelerated. Here we use axisymmetric relativistic MHD simulations to obtain constraints on the acceleration site(s) of particles of different energies in the Crab Nebula. Various scenarios are considered for the injection of particles responsible for synchrotron radiation in the different frequency bands, radio, optical and X-rays. The resulting emission properties are compared with available data on the multi wavelength time variability of the inner nebula. Our main result is that the X-ray emitting particles are accelerated in the equatorial region of the pulsar wind. Possible implications on the nature of the acceleration mechanism are discussed.

Magnetic fields play a critical role in the phenomenology of neutron stars. There is virtually no observable aspect which is not governed by them. Despite this, only recently efforts have been done to model magnetic fields in the correct general relativistic regime, characteristic of these compact objects. In this work we present, for the first time a comprehensive and detailed parameter study, in general relativity, of the role that the current distribution, and the related magnetic field structure, have in determining the precise structure of neutron stars. In particular, we show how the presence of localized currents can modify the field strength at the stellar surface, and we look for general trends, both in terms of energetic properties, and magnetic field configurations. Here we verify that, among other things, for a large class of different current distributions the resulting magnetic configurations are always dominated by the poloidal component of the current.

Soft Gamma-Ray Repeaters and Anomalous X-Ray Pulsars are extreme manifestations of the most magnetized neutron stars: magnetars. The phenomenology of their emission and spectral properties strongly support the idea that the magnetospheres of these astrophysical objects are tightly twisted in the vicinity of the star. Previous studies on equilibrium configurations have so far focused on either the internal or the external magnetic field configuration, without considering a real coupling between the two fields. Here we investigate numerical equilibrium models of magnetized neutron stars endowed with a confined twisted magnetosphere, solving the general relativistic Grad-Shafranov equation both in the interior and in the exterior of the compact object. A comprehensive study of the parameters space is provided to investigate the effects of different current distributions on the overall magnetic field structure.

The long-term evolution of large-amplitude Alfven waves propagating in the solar wind is investigated by performing two-dimensional MHD simulations within the expanding box model. The linear and nonlinear phases of the parametric decay instability are studied for both circularly polarized waves in parallel propagation and for arc-polarized waves in oblique propagation. The non-monochromatic case is also considered. In the oblique case, the direct excitation of daughter modes transverse to the local background field is found for the first time in an expanding environment, and this transverse cascade seems to be favored for monochromatic mother waves. The expansion effect reduces the instability growth rate, and it can even suppress its onset for the lowest frequency modes considered here, possibly explaining the persistence of these outgoing waves in the solar wind.

High angular resolution spectra obtained with the Hubble Space Telescope Imaging Spectrograph (HST/STIS) provide rich morphological and kinematical information about the stellar jet phenomenon, which allows us to test theoretical models efficiently. In this work, numerical simulations of stellar jets in the propagation region are executed with the PLUTO code, by adopting inflow conditions that arise from former numerical simulations of magnetized outflows, accelerated by the disk-wind mechanism in the launching region. By matching the two regions, information about the magneto-centrifugal accelerating mechanism underlying a given astrophysical object can be extrapolated by comparing synthetic and observed position-velocity diagrams (PVDs). We show that quite different jets, like those from the young T Tauri stars DG-Tau and RW-Aur, may originate from the same disk-wind model for different configurations of the magnetic field at the disk surface. This result supports the idea that all the observed jets may be generated by the same mechanism.

Extremely magnetized neutron stars with magnetic fields as strong as $\sim 10^{15-16}$ G, or magnetars, have received considerable attention in the last decade due to their identification as a plausible source for Soft Gamma Repeaters and Anomalous X-ray Pulsars. Moreover, this class of compact objects has been proposed as a possible engine capable of powering both Long and Short Gamma-Ray Bursts, if the rotation period in their formation stage is short enough (~1 ms). Such strong fields are expected to induce substantial deformations of the star and thus to produce the emission of gravitational waves. Here we investigate, by means of numerical modeling, axisymmetric static equilibria of polytropic and strongly magnetized stars in full general relativity, within the ideal magneto-hydrodynamic regime. The eXtended Conformally Flat Condition (XCFC) for the metric is assumed, allowing us to employ the techniques introduced for the X-ECHO code [Bucciantini & Del Zanna, 2011, Astron. Astrophys. 528, A101], proven to be accurate, efficient, and stable. The updated XNS code for magnetized neutron star equilibria is made publicly available for the community (see www.arcetri.astro.it/science/ahead/XNS). Several sequences of models are here retrieved, from the purely toroidal (resolving a controversy in the literature) or poloidal cases, to the so-called twisted torus mixed configurations, expected to be dynamically stable, which are solved for the first time in the non-perturbative regime.

We present the first kinematic study of an $\alpha\Omega$-dynamo in the General Relativistic Magneto-HydroDynamics (GRMHD) regime, applied to thick disks orbiting around Kerr black holes and using a fully covariant mean field dynamo closure for the Ohm law. We show that the $\alpha\Omega$-dynamo mechanism leads to a continuous exponential growth of the magnetic field within the disk and to the formation of dynamo waves drifting away or toward the equatorial plane. Since the evolution of the magnetic field occurs qualitatively in the same fashion as in the Sun, we present also butterfly diagrams that characterize our models and show the establishment of an additional timescale, which depends on the microscopic properties of the turbulent motions, possibly providing an alternative explanation to periodicities observed in many high-energy astrophysical sources where accretion onto a rotating black hole is believed to operate.

High-energy phenomena in astrophysics involve quite generally a combination of relativistic motions and strong gravity. The simultaneous solution of Einstein equations and General Relativistic MHD equations is thus necessary to model with accuracy such phenomena. The so-called Conformally Flat Condition (CFC) allows a simplified treatment of Einstein equations, that can be particularly efficient in those contexts where gravitational wave emission is negligible, like core-collapse, or the formation/evolution of neutron stars. We have developed a set of codes to model axisymmetric MHD flows, in General Relativity, where the solution of Einstein equations is achieved with a semi-spectral scheme. Here, we will show how this framework is particularly well suited to investigate neutron star equilibrium models in the presence of strong magnetic fields and we will present the XNS code, that has been recently developed and here updated to treat poloidal and mixed configurations.

We present the first astrophysical application of the ECHO code in its recent version supplemented by a generalized Ohm law, namely a kinematic study of dynamo effects in thick accretion disks. High-order \emph{implicit-explicit} Runge-Kutta time-stepping routines are implemented and validated within $3+1$ \emph{General Relativistic MagnetoHydroDynamics} (GRMHD). The scheme is applied to a differentially rotating torus orbiting a Kerr black hole, where the mean-field dynamo process leads to strong amplification of seed magnetic fields. We show that the interplay between the toroidal and poloidal components occurs qualitatively in the same fashion as in the Sun, butterfly diagrams are reproduced, and a typical time-scale for the field evolution is found, depending on the dynamo and resistivity numbers, which could explain periodicities as observed in several accreting systems.

Magnetic fields play a crucial role in many astrophysical scenarios and, in particular, are of paramount importance in the emission mechanism and evolution of Neutron Stars (NSs). To understand the role of the magnetic field in compact objects it is important to obtain, as a first step, accurate equilibrium models for magnetized NSs. Using the conformally flat approximation we solve the Einstein's equations together with the GRMHD equations in the case of a static axisymmetryc NS taking into account different types of magnetic configuration. This allows us to investigate the effect of the magnetic field on global properties of NSs such as their deformation.

In recent years, it has become a well-established paradigm that many aspects of the physics of Pulsar Wind Nebulae (PWNe) can be fully accounted for within a relativistic MHD description. Numerical simulations have proven extremely successful in reproducing the X-ray morphology of the Crab Nebula, down to very fine detail. Radio emission, instead, is currently one of the most obscure aspects of the physics of these objects, and one that holds important information about pulsar properties and their role as antimatter factories. Here we address the question of radio emission morphology and integrated spectrum from the Crab Nebula, by using for the first time an axisymmetric dynamical model with parameters chosen to best reproduce its X-ray morphology. Based on our findings we discuss constraints on the origin of the radio emitting particles.

We present ECHO-QGP, a numerical code for $(3+1)$-dimensional relativistic viscous hydrodynamics designed for the modeling of the space-time evolution of the matter created in high energy nuclear collisions. The code has been built on top of the \emph{Eulerian Conservative High-Order} astrophysical code for general relativistic magneto-hydrodynamics [\emph{Del Zanna et al., Astron. Astrophys. 473, 11, 2007}] and here it has been upgraded to handle the physics of the Quark-Gluon Plasma. ECHO-QGP features second-order treatment of causal relativistic viscosity effects in both Minkowskian or Bjorken coordinates; partial or complete chemical equilibrium of hadronic species before kinetic freeze-out; initial conditions based on the optical Glauber model, including a Monte-Carlo routine for event-by-event fluctuating initial conditions; a freeze-out procedure based on the Cooper-Frye prescription. The code is extensively validated against several test problems and results always appear accurate, as guaranteed by the combination of the conservative (shock-capturing) approach and the high-order methods employed. ECHO-QGP can be extended to include evolution of the electromagnetic fields coupled to the plasma.

The powerful high-energy phenomena typically encountered in astrophysics invariably involve physical engines, like neutron stars and black hole accretion disks, characterized by a combination of highly magnetized plasmas, strong gravitational fields, and relativistic motions. In recent years numerical schemes for General Relativistic MHD (GRMHD) have been developed to model the multidimensional dynamics of such systems, including the possibility of an evolving spacetime. Such schemes have been also extended beyond the ideal limit including the effects of resistivity, in an attempt to model dissipative physical processes acting on small scales (sub-grid effects) over the global dynamics. Along the same lines, magnetic fields could be amplified by the presence of turbulent dynamo processes, as often invoked to explain the high values of magnetization required in accretion disks and neutron stars. Here we present, for the first time, a further extension to include the possibility of a mean-field dynamo action within the framework of numerical 3+1 (resistive) GRMHD. A fully covariant dynamo closure is proposed, in analogy with the classical theory, assuming a simple alpha-effect in the comoving frame. Its implementation into a finite-difference scheme for GRMHD in dynamical spacetimes [the X-ECHO code: (Bucciantini and Del Zanna 2011)] is described, and a set of numerical test is presented and compared with analytical solutions wherever possible.

We present a new numerical code, X-ECHO, for general relativistic magnetohydrodynamics (GRMHD) in dynamical spacetimes. This is aimed at studying astrophysical situations where strong gravity and magnetic fields are both supposed to play an important role, such as for the evolution of magnetized neutron stars or for the gravitational collapse of the magnetized rotating cores of massive stars, which is the astrophysical scenario believed to eventually lead to (long) GRB events. The code is based on the extension of the Eulerian conservative high-order (ECHO) scheme [Del Zanna et al., A&A 473, 11 (2007)] for GRMHD, here coupled to a novel solver for the Einstein equations in the extended conformally flat condition (XCFC). We fully exploit the 3+1 Eulerian formalism, so that all the equations are written in terms of familiar 3D vectors and tensors alone, we adopt spherical coordinates for the conformal background metric, and we consider axisymmetric spacetimes and fluid configurations. The GRMHD conservation laws are solved by means of shock-capturing methods within a finite-difference discretization, whereas, on the same numerical grid, the Einstein elliptic equations are treated by resorting to spherical harmonics decomposition and solved, for each harmonic, by inverting band diagonal matrices. As a side product, we build and make available to the community a code to produce GRMHD axisymmetric equilibria for polytropic relativistic stars in the presence of differential rotation and a purely toroidal magnetic field. This uses the same XCFC metric solver of the main code and has been named XNS. Both XNS and the full X-ECHO codes are validated through several tests of astrophysical interest.

We investigate the dynamics of a circumbinary disc that responds to the loss of mass and to the recoil velocity of the black hole produced by the merger of a binary system of supermassive black holes. We perform the first two-dimensional general relativistic hydrodynamics simulations of \textit{extended} non-Keplerian discs and employ a new technique to construct a "shock detector", thus determining the precise location of the shocks produced in the accreting disc by the recoiling black hole. In this way we can study how the properties of the system, such as the spin, mass and recoil velocity of the black hole, affect the mass accretion rate and are imprinted on the electromagnetic emission from these sources. We argue that the estimates of the bremsstrahlung luminosity computed without properly taking into account the radiation transfer yield cooling times that are unrealistically short. At the same time we show, through an approximation based on the relativistic isothermal evolution, that the luminosity produced can reach a peak value above $L \simeq 10^{43} \ {\rm erg/s} $ at about $\sim 30\,{\rm d}$ after the merger of a binary with total mass $M\simeq 10^6 M_\odot$ and persist for several days at values which are a factor of a few smaller. If confirmed by more sophisticated calculations such a signal could indeed lead to an electromagnetic counterpart of the merger of binary black-hole system.

The main goal of our present work is to provide, for the first time, a simple computational tool that can be used to compute the brightness, the spectral index, the polarization, the time variability and the spectrum of the non-thermal light (both synchrotron and inverse Compton, IC) associated with the plasma dynamics resulting from given relativistic magnetohydrodynamics (RMHD) simulations. The proposed method is quite general, and can be applied to any scheme for RMHD and to all non-thermal emitting sources, e.g. pulsar wind nebulae (PWNe), and in particular to the Crab Nebula (CN) as in the present proceeding. Here only the linear optical and X-ray polarization that characterizes the PWNe synchrotron emission is analyzed in order to infer information on the inner bulk flow structure, to provide a direct investigation of the magnetic field configuration, in particular the presence and the strength of a poloidal component, and to understand the origin of some emitting features, such as the knot, whose origins are still uncertain. The inverse Compton radiation is examined to disentangle the different contributions to radiation from the magnetic field and the particle energy distribution function, and to search for a possible hadronic component in the emitting PWN, and thus for the presence of ions in the wind. If hadronic radiation was found in a PWN, young supernova remnants would provide a natural birth-place of the cosmic-rays (CRs) up to the so-called knee in the CR spectrum.

We use ideal axisymmetric relativistic magnetohydrodynamic simulations to calculate the spindown of a newly formed millisecond, B ~ 10^{15} G, magnetar and its interaction with the surrounding stellar envelope during a core-collapse supernova (SN) explosion. The mass, angular momentum, and rotational energy lost by the neutron star are determined self-consistently given the thermal properties of the cooling neutron star's atmosphere and the wind's interaction with the surrounding star. The magnetar drives a relativistic magnetized wind into a cavity created by the outgoing SN shock. For high spindown powers (~ 10^{51}-10^{52} ergs/s), the magnetar wind is super-fast at almost all latitudes, while for lower spindown powers (~ 10^{50} erg/s), the wind is sub-fast but still super-Alfvenic. In all cases, the rates at which the neutron star loses mass, angular momentum, and energy are very similar to the corresponding free wind values (<~ 30% differences), in spite of the causal contact between the neutron star and the stellar envelope. In addition, in all cases that we consider, the magnetar drives a collimated (~5-10 deg.) relativistic jet out along the rotation axis of the star. Nearly all of the spindown power of the neutron star escapes via this polar jet, rather than being transferred to the more spherical SN explosion. The properties of this relativistic jet and its expected late-time evolution in the magnetar model are broadly consistent with observations of long duration gamma-ray bursts (GRBs) and their associated broad-lined Type Ic SN.

We study the evolution of an accelerating hyperrelativistic shock under the presence of upstream inhomogeneities wrinkling the discontinuity surface. The investigation is conducted by means of numerical simulations using the PLUTO code for astrophysical fluid dynamics. The reliability and robustness of the code are demonstrated against well known results coming from the linear perturbation theory. We then follow the nonlinear evolution of two classes of perturbing upstream atmospheres and conclude that no lasting wrinkle can be preserved indefinitely by the flow. Finally we derive analytically a description of the geometrical effects of a turbulent upstream ambient on the discontinuity surface.

In this paper we complete the set of diagnostic tools for synchrotron emitting sources presented by Del Zanna et al. (Astron. Astrophys. 453, 621, 2006) with the computation of inverse Compton radiation from the same relativistic particles. Moreover we investigate, for the first time, the gamma-ray emission properties of Pulsar Wind Nebulae in the light of the axisymmetric jet-torus scenario. The method consists in evolving the relativistic MHD equations and the maximum energy of the emitting particles. The particle energy distribution function is split in two components: the radio one connected to a relic population born at the outburst of the supernova and the other associated to the wind population continuously accelerated at the termination shock and emitting up to the gamma-ray band. We consider the general Klein-Nishina cross section and three different photon targets: the nebular synchrotron photons, far-infrared thermal ones and the cosmic microwave background. The overall synchrotron spectrum is fitted assuming an excess of injected particles and a steeper power law with respect to previous models. The TeV emission has the correct shape but is in excess of the data. This is due to the nebular magnetic field structure as obtained by the simulations. The jet-torus morphology is visible in high-resolution gamma-ray synthetic maps too. We present a preliminary exploration of time variability in the X and gamma-ray bands.

We present a complete set of diagnostic tools aimed at reproducing synthetic non-thermal (synchrotron and/or Inverse Compton, IC) emissivity, integrated flux energy, polarization and spectral index simulated maps in comparison to observations. The time dependent relativistic magnetohydrodynamic (RMHD) equations are solved with a shock capturing code together with the evolution of the maximum particles energy. Applications to Pulsar Wind Nebulae (PWNe) are shown.

Elongated jets from young stellar objects typically present a nodular structure, formed by a chain of bright knots of enhanced emission with individual proper motions. Though it is generally accepted that internal shocks play an important role in the formation and dynamics of such structures, their precise origin and the mechanisms behind the observed proper motions is still a matter of debate. Our goal is to study numerically the origin, dynamics, and emission properties of such knots. Axisymmetric simulations are performed with a shock-capturing code for gas dynamics, allowing for molecular, atomic, and ionized hydrogen in non-equilibrium concentrations subject to ionization/recombination processes. Radiative losses in SII lines are computed, and the resulting synthetic emission maps are compared with observations. We show that a pattern of regularly spaced internal oblique shocks, characterized by individual proper motions, is generated by the pressure gradient between the propagating jet and the time variable external cocoon. In the case of under-expanded, light jets the resulting emission knots are found to move downstream with the jet flow, with increasing velocity and decaying brightness toward the leading bow shock. This suggests that the basic properties of the knots observed in stellar jets can be reproduced even without invoking ad hoc pulsating conditions at the jet inlet, though an interplay between the two scenarios is certainly possible.

We present a new numerical code, ECHO, based on an Eulerian Conservative High Order scheme for time dependent three-dimensional general relativistic magnetohydrodynamics (GRMHD) and magnetodynamics (GRMD). ECHO is aimed at providing a shock-capturing conservative method able to work at an arbitrary level of formal accuracy (for smooth flows), where the other existing GRMHD and GRMD schemes yield an overall second order at most. Moreover, our goal is to present a general framework, based on the 3+1 Eulerian formalism, allowing for different sets of equations, different algorithms, and working in a generic space-time metric, so that ECHO may be easily coupled to any solver for Einstein's equations. Various high order reconstruction methods are implemented and a two-wave approximate Riemann solver is used. The induction equation is treated by adopting the Upwind Constrained Transport (UCT) procedures, appropriate to preserve the divergence-free condition of the magnetic field in shock-capturing methods. The limiting case of magnetodynamics (also known as force-free degenerate electrodynamics) is implemented by simply replacing the fluid velocity with the electromagnetic drift velocity and by neglecting the matter contribution to the stress tensor. ECHO is particularly accurate, efficient, versatile, and robust. It has been tested against several astrophysical applications, including a novel test on the propagation of large amplitude circularly polarized Alfven waves. In particular, we show that reconstruction based on a Monotonicity Preserving filter applied to a fixed 5-point stencil gives highly accurate results for smooth solutions, both in flat and curved metric (up to the nominal fifth order), while at the same time providing sharp profiles in tests involving discontinuities.

We present here a series of numerical simulations of the development of Kelvin-Helmholtz instability in a relativistically hot plasma. The physical parameters in the unperturbed state are chosen to be representative of local conditions encountered in Pulsar Wind Nebulae (PWNe), with a main magnetic field perpendicular to a mildly relativistic shear layers. By using a numerical code for Relativistic MHD, we investigate the effect of an additional magnetic field component aligned with the shear velocity, and we follow the evolution of the instability to the saturation and turbulent regimes. Based on the resulting flow structure, we then compute synchrotron maps in order to evaluate the signature of Kelvin-Helmholtz instability on the emission and we investigate how the time scale and the amplitude of the synchrotron modulations depend on shear velocity and magnetic field. Finally we compare our results to the observed variable features in the Crab Nebula. We show that the Kelvin-Helmholtz instability cannot account for the wisps variability, but it might be responsible for the time dependent filamentary structure observed in the main torus.

A complete set of diagnostic tools aimed at producing synthetic synchrotron emissivity, polarization, and spectral index maps from relativistic MHD simulations is presented. As a first application we consider here the case of the emission from Pulsar Wind Nebulae (PWNe). The proposed method is based on the addition, on top of the basic set of MHD equations, of an extra equation describing the evolution of the maximum energy of the emitting particles. This equation takes into account adiabatic and synchrotron losses along streamlines for the distribution of emitting particles and its formulation is such that it is easily implemented in any numerical scheme for relativistic MHD. Application to the axisymmetric simulations of PWNe, analogous to those described by Del Zanna et al. (2004, A&A, 421, 1063), allows direct comparison between the numerical results and observations of the inner structure of the Crab Nebula, and similar objects, in the optical and X-ray bands. We are able to match most of the observed features typical of PWNe, like the equatorial torus and the polar jets, with velocities in the correct range, as well as finer emission details, like arcs, rings and the bright knot, that turn out to arise mainly from Doppler boosting effects. Spectral properties appear to be well reproduced too: detailed spectral index maps are produced for the first time and show softening towards the PWN outer borders, whereas spectral breaks appear in integrated spectra. The emission details are found to strongly depend on both the average wind magnetization (here approximately 2%), and on the magnetic field shape.

We solve for the time-dependent dynamics of axisymmetric, general relativistic magnetohydrodynamic winds from rotating neutron stars. The mass loss rate is obtained self-consistently as a solution to the MHD equations, subject to a finite thermal pressure at the stellar surface. We consider both monopole and dipole magnetic field geometries and we explore the parameter regime extending from low magnetization (low-sigma_o), almost thermally-driven winds to high magnetization (high-sigma_o), relativistic Poynting-flux dominated outflows. We compute the angular momentum and rotational energy loss rates as a function of sigma_o and compare with analytic expectations from the classical theory of pulsars and magnetized stellar winds. In the case of the monopole, our high-sigma_o calculations asymptotically approach the analytic force-free limit. If we define the spindown rate in terms of the open magnetic flux, we similarly reproduce the spindown rate from recent force-free calculations of the aligned dipole. However, even for sigma_o as high as ~20, we find that the location of the Y-type point (r_Y), which specifies the radius of the last closed field line in the equatorial plane, is not the radius of the light cylinder R_L = c/omega (R = cylindrical radius), as has previously been assumed in most estimates and force-free calculations. Instead, although the Alfven radius at intermediate latitudes quickly approaches R_L as sigma_o exceeds unity, r_Y remains significantly less than R_L. Because r_Y < R_L, our calculated spindown rates thus exceed the classic ``vacuum dipole'' rate. We discussthe implications of our results for models of rotation-powered pulsars and magnetars, both in their observed states and in their hypothesized rapidly rotating initial state.

Conservative numerical schemes for general relativistic magnetohydrodynamics (GRMHD) require a method for transforming between ``conserved'' variables such as momentum and energy density and ``primitive'' variables such as rest-mass density, internal energy, and components of the four-velocity. The forward transformation (primitive to conserved) has a closed-form solution, but the inverse transformation (conserved to primitive) requires the solution of a set of five nonlinear equations. Here we discuss the mathematical properties of the inverse transformation and present six numerical methods for performing the inversion. The first method solves the full set of five nonlinear equations directly using a Newton-Raphson scheme and a guess from the previous timestep. The other methods reduce the five nonlinear equations to either one or two nonlinear equations that are solved numerically. Comparisons between the methods are made using a survey over phase space, a two-dimensional explosion problem, and a general relativistic MHD accretion disk simulation. The run-time of the methods is also examined. Code implementing the schemes is available for download on the web.

We present here the first effort to compute synthetic synchrotron polarization maps of Pulsar Wind Nebulae (PWNe). Our goal is to highlight how polarization can be used as an additional diagnostic tool for the flow structure in the inner regions of these nebulae. Recent numerical simulations suggest the presence of flow velocities ~0.5 c in the surroundings of the termination shock, where most of the high energy emission comes from. We construct polarization maps taking into account relativistic effects like Doppler boosting and position angle swing. The effect of different bulk velocities is clarified with the help of a toy-model consisting of a uniformly emitting torus. We also present a map based on recent numerical simulations of the entire nebula and compare it with presently available data. The comparison with upcoming high resolution observations could provide new insight into the inner structure of the nebula and put constraints on the geometrical properties of the magnetic field.

Pulsar bow-shock nebulae are a class of pulsar wind nebulae (PWNe) that form when the pulsar wind is confined by the ram pressure of the ambient medium, and are usually associated with old pulsars, that have already emerged from the progenitor Supernova Remnant (SNR). Until a few years ago these nebulae were mainly observed as \halpha sources; recently, also non-thermal emission has been detected. This is the signature of accelerated particles gyrating in a magnetic field. In the same way as \halpha radiation is a tool for studying the layer of shocked Interstellar Medium (ISM), the non-thermal radiation might be used to infer the properties of the shocked pulsar wind. However theoretical and numerical models have been presented so far only in the hydrodynamical (HD) regime, while in order to properly model the internal flow structure and the emission properties of these nebulae a magnetohydrodynamical (MHD) treatment is required. We present here relativistic MHD (RMHD) axisymmetric simulations of pulsar wind bow-shock nebulae. The structure and fluid dynamics of such objects is investigated for various values of the pulsar wind magnetization. Simulated synchrotron maps are computed and comparison of the emission pattern with observations is discussed.

The propagation and the evolution of Alfvenic pulses in the solar coronal arcades is investigated by means of MHD numerical simulations. Significant transverse oscillations in coronal loops, triggered by nearby flare events, are often measured in EUV lines and are generally interpreted as standing kink modes. However, the damping times of these oscillations are typically very short (from one to a few periods) and the physical mechanism responsible for the decay is still a matter of debate. Moreover, the majority of the observed cases actually appears to be better modeled by propagating, rather than standing, modes. Here we perform 2.5-D compressible MHD simulations of impulsively generated Alfven waves propagating in a potential magnetic arcade (assumed as a simplified 2-D loop model), taking into account the stratification of the solar atmosphere with height from the photosphere to the corona. The results show a strong spreading of the initially localized pulses along the loop, due to the variations in the Alfven velocity with height, and correspondingly an efficient damping of the amplitude of the oscillations. We believe that simple explanations based on the effects of wave propagation in highly inhomogeneous media may apply to the majority of the reported cases, and that variations of the background density and Alfven speed along the loop should be considered as key ingredients in future models.

We present a numerical investigation of the development of Rayleigh-Taylor instability at the interface between an expanding Pulsar Wind Nebula and its surrounding Supernova Remnant. These systems have long been thought to be naturally subject to this kind of instability, given their expansion behavior and the density jump at the contact discontinuity. High resolution images of the Crab Nebula at optical frequencies show the presence of a complex network of line-emitting filaments protruding inside the synchrotron nebula. These structures are interpreted as the observational evidence that Rayleigh-Taylor instability is in fact at work. The development of this instability in the regime appropriate to describe Supernova Remnant-Pulsar Wind Nebula systems is non-trivial. The conditions at the interface are likely close to the stability threshold, and the inclusion of the nebular magnetic field, which might play an important role in stabilizing the system, is essential to the modeling. If Rayleigh-Taylor features can grow efficiently a mixing layer in the outer portion of the nebula might form where most of the supernova material is confined. When a magnetic field close to equipartition is included we find that the interface is stable, and that even a weaker magnetic field affects substantially the growth and shape of the fingers.

The structure and the evolution of Pulsar Wind Nebulae (PWNe) are studied by means of two-dimensional axisymmetric relativistic magnetohydrodynamic (RMHD) simulations. After the first imaging of the Crab Nebula with Chandra, a growing number of objects has been found to show in the X-rays spatial features such as rings and jets, that clearly cannot be accounted for within the standard framework of one-dimensional semi-analytical models. The most promising explanation suggested so far is based on the combined effects of the latitude dependence of the pulsar wind energy flux, shaping the wind termination shock and naturally providing a higher equatorial emission, and of the wind magnetization, likely responsible for the jet collimation by hoop stresses downstream of the shock. This scenario is investigated here by following the evolution of a PWN interacting with the confining Supernova Remnant (SNR), from the free expansion to the beginning of the reverberation phase. Our results confirm the oblate shape of the wind termination shock and the formation of a polar jet with supersonic velocities (v~0.5-0.7 c) for high enough values of the equatorial wind magnetization parameter (sigma~0.01).

We present high resolution spherically symmetric relativistic magnetohydrodynamical simulations of the evolution of a pulsar wind nebula inside the free expanding ejecta of the supernova progenitor. The evolution is followed starting from a few years after the supernova explosion and up to an age of the remnant of 1500 years. We consider different values of the pulsar wind magnetization parameter and also different braking indices for the spin-down process. We compare the numerical results with those derived through an approximate semi-analytical approach that allows us to trace the time evolution of the positions of both the pulsar wind termination shock and the contact discontinuity between the nebula and the supernova ejecta. We also discuss, whenever a comparison is possible, to what extent our numerical results agree with former self-similar models, and how these models could be adapted to take into account the temporal evolution of the system. The inferred magnetization of the pulsar wind could be an order of magnitude lower than that derived from time independent analytic models.

We present a general framework to design Godunov-type schemes for multidimensional ideal magnetohydrodynamic (MHD) systems, having the divergence-free relation and the related properties of the magnetic field B as built-in conditions. Our approach mostly relies on the 'Constrained Transport' (CT) discretization technique for the magnetic field components, originally developed for the linear induction equation, which assures div(B)=0 and its preservation in time to within machine accuracy in a finite-volume setting. We show that the CT formalism, when fully exploited, can be used as a general guideline to design the reconstruction procedures of the B vector field, to adapt standard upwind procedures for the momentum and energy equations, avoiding the onset of numerical monopoles of O(1) size, and to formulate approximate Riemann solvers for the induction equation. This general framework will be named here 'Upwind Constrained Transport' (UCT). To demonstrate the versatility of our method, we apply it to a variety of schemes, which are finally validated numerically and compared: a novel implementation for the MHD case of the second order Roe-type positive scheme by Liu and Lax (J. Comp. Fluid Dynam. 5, 133, 1996), and both the second and third order versions of a central-type MHD scheme presented by Londrillo and Del Zanna (Astrophys. J. 530, 508, 2000), where the basic UCT strategies have been first outlined.

Pulsars, formed during supernova explosions, are known to be sources of relativistic magnetized winds whose interaction with the expanding supernova remnants (SNRs) gives rise to a pulsar wind nebula (PWN). We present spherically symmetric relativistic magnetohydrodynamics (RMHD) simulations of the interaction of a pulsar wind with the surrounding SNR, both in particle and magnetically dominated regimes. As shown by previous simulations, the evolution can be divided in three phases: free expansion, a transient phase characterized by the compression and reverberation of the reverse shock, and a final Sedov expansion. The evolution of the contact discontinuity between the PWN and the SNR (and consequently of the SNR itself) is almost independent of the magnetization of the nebula as long as the total (magnetic plus particle) energy is the same. However, a different behaviour of the PWN internal structure is observable during the compression-reverberation phase, depending on the degree of magnetization=2E The simulations were performed using the third order conservative scheme by Del Zanna et al. (2003).

A third order shock-capturing numerical scheme for three-dimensional special relativistic magnetohydrodynamics (3-D RMHD) is presented and validated against several numerical tests. The simple and efficient central scheme described in Paper I (Del Zanna and Bucciantini, Astron. Astrophys., 390, 1177--1186, 2002) for relativistic hydrodynamics is here extended to the magnetic case by following the strategies prescribed for classical MHD by Londrillo and Del Zanna (Astrophys. J., 530, 508--524, 2000). The scheme avoids completely spectral decomposition into characteristic waves, computationally expensive and subject to many degenerate cases in the magnetic case, while it makes use of a two-speed Riemann solver that just require the knowledge of the two local fast magnetosonic velocities. Moreover, the onset of spurious magnetic monopoles, which is a typical problem for multi-dimensional MHD upwind codes, is prevented by properly taking into account the solenoidal constraint and the specific antisymmetric nature of the induction equation. Finally, the extension to generalized orthogonal curvilinear coordinate systems is included, thus the scheme is ready to incorporate general relativistic (GRMHD) effects.

Multidimensional shock-capturing numerical schemes for special relativistic hydrodynamics (RHD) are computationally more expensive than their correspondent Euler versions, due to the nonlinear relations between conservative and primitive variables and to the consequent complexity of the Jacobian matrices (needed for the spectral decomposition in most of the approximate Riemann solvers of common use). Here an efficient and easy-to-implement three-dimensional (3-D) shock-capturing scheme for ideal RHD is presented. Based on the algorithms developed by P. Londrillo and L. Del Zanna ({\em Astrophys. J.} 530, 508-524, 2000) for the non-relativistic magnetohydrodynamic (MHD) case, and having in mind its relativistic MHD extension (to appear in a forthcoming paper), the scheme uses high order (third) Convex Essentially Non-Oscillatory (CENO) finite difference interpolation routines and central-type averaged Riemann solvers, which do not make use of time-consuming characteristic decomposition. The scheme is very efficient and robust, and it gives results comparable to those obtained with more sophisticated algorithms, even in ultrarelativistic multidimensional test problems.

A general method for constructing high order upwind schemes for multidimensional magnetohydrodynamics (MHD), having as a main built-in condition the divergence-free constraint $\divb=0$ for the magnetic field vector $\bb$, is proposed. The suggested procedure is based on {\em consistency} arguments, by taking into account the specific operator structure of MHD equations with respect to the reference Euler equations of gas-dynamics. This approach leads in a natural way to a staggered representation of the $\bb$ field numerical data where the divergence-free condition in the cell-averaged form, corresponding to second order accurate numerical derivatives, is exactly fulfilled. To extend this property to higher order schemes, we then give general prescriptions to satisfy a $(r+1)^{th}$ order accurate $\divb=0$ relation for any numerical $\bb$ field having a $r^{th}$ order interpolation accuracy. Consistency arguments lead also to a proper formulation of the upwind procedures needed to integrate the induction equations, assuring the exact conservation in time of the divergence-free condition and the related continuity properties for the $\bb$ vector components. As an application, a third order code to simulate multidimensional MHD flows of astrophysical interest is developed using ENO-based reconstruction algorithms. Several test problems to illustrate and validate the proposed approach are finally presented.