Dust radiative transfer modelling of the infrared ring around the magnetar SGR 1900$+$14
G. Natale N. Rea D. Lazzati R. Perna D. F. Torres J. M. Girart




G. Natale
gnatale@uclan.ac.uk


Jeremiah Horrocks Institute, University of Central Lancashire, Preston, PR1 2HE, UK


Institute of Space Sciences (IEEC–CSIC), Campus UAB, Carrer de Can Magrans S/N, 08193 Barcelona, Spain.
Anton Pannekoek Institute for Astronomy, University of Amsterdam, Postbus 94249, NL–1090 GE Amsterdam, the Netherlands.


Department of Physics, Oregon State University, 301 Weniger Hall, Corvallis, OR 97331, USA


Department of Physics and Astronomy, Stony Brook University, Stony Brook, NY, 11794, USA


Institute of Space Sciences (IEEC–CSIC), Campus UAB, Carrer de Can Magrans S/N, 08193 Barcelona, Spain.
Institució Catalana de Recerca i Estudis Avançats (ICREA), E-08010 Barcelona, Spain


Institute of Space Sciences (IEEC–CSIC), Campus UAB, Carrer de Can Magrans S/N, 08193 Barcelona, Spain.




A peculiar infrared ring-like structure was discovered by  Spitzer  around the strongly magnetised neutron star  SGR 1900+14. This infrared structure was suggested to be due to a dust-free cavity, produced by the SGR Giant Flare occurred in 1998, and kept illuminated by surrounding stars. Using a 3D dust radiative transfer code, we aimed at reproducing the emission morphology and the integrated emission flux of this structure assuming different spatial distributions and densities for the dust, and different positions for the illuminating stars.  We found that a dust-free ellipsoidal cavity can reproduce the shape, flux, and spectrum of the ring-like infrared emission, provided that the illuminating stars are inside the cavity and that the interstellar medium has high gas density (n_H∼1000 cm^-3). We further constrain the emitting region to have a sharp inner boundary and to be significantly extended in the radial direction, possibly even just a cavity in a smooth molecular cloud. We discuss possible scenarios for the formation of the dustless cavity and the particular geometry that allows it to be IR-bright.












§ INTRODUCTION
 


Strongly magnetized neutron stars <cit.> are extremely powerful X-ray and soft gamma-ray emitters, in particular under the form of large flares. These flares might reach luminosities that in our Galaxy second only Supernova explosions. In particular, magnetars emit a large variety of flares and outbursts on timescales from fraction of seconds to years, on a vast range of luminosities from 10^38 to ∼10^47<cit.>. The most energetic events they ever emitted, called Giant Flares, have been detected three times in the past few decades, from three magnetars: the Soft Gamma-ray Repeaters (SGR) 0526-66 on 1979 March 5 <cit.>, SGR 1900+14 on 1998 August 27 <cit.> and the last and more energetic one on 2004 December 27 from SGR 1806-20 <cit.>. All of them are characterized by a very luminous initial spike (∼10^44-47), lasting less than a second, which decays rapidly into a softer tail (modulated at the neutron star spin period) lasting several hundreds of seconds (with luminosities of ∼10^43). The nature of the steady and flaring high energy emission from these sources has been intriguing all along. In fact, magnetar's X-ray luminosity is in general too high to be produced by the pulsar rotational energy losses alone, as for more common isolated radio pulsars, and the lack of any companion star excludes an accretion scenario. It is now well established that the peculiarities of these extreme highly magnetized objects (10^14-15 Gauss) are related to the strength and instability of their magnetic field, that at times might stress the stiff neutron star crust <cit.>, might rearrange itself locally in small twisted bundles, or disrupt and reconnect higher up in the magnetosphere producing large ejections of particles. The ages of the ∼25 magnetars known <cit.>, derived from their rotational properties (t_c∼Ṗ/P), indicate a young population, typically a few thousand years old. In three or four cases there are reasonably accepted associations with Supernova Remnants <cit.>, as well as massive star clusters <cit.>.  





  

 
SGR1900+14 is one of the youngest magnetars known (∼1 kyr).  A very prolific burster <cit.>, it is embedded in a cluster of very massive stars <cit.> and it is one of the three magnetars which has shown a Giant Flare. This magnetar was observed using all three instruments onboard the NASA Spitzer Space Telescope in 2005 and 2007 <cit.>.
Surprisingly, these observations have revealed a prominent ring-like structure (see Figure 1) in the 16μm and 24μm wave-bands, not detected in the 3.6–8.0 μm observations. A formal elliptical fit to the ring indicates semi-major and semi-minor axes of angular lengths ∼36” and ∼19”, respectively, centred at the position of the magnetar SGR1900+14. No equivalent feature was observed at radio or X-ray wavelengths (L_ 332 MHz≤ 2.7 × 10^29 d_12.5^2  erg s^-1, and L_ 2-10 keV≤ 1.8 × 10^33 d_12.5^2  erg s^-1; with d_12.5 being the distance in units of 12.5 kpc; <cit.>).

The Spitzer images are dominated by the bright emission from two nearby M supergiants that mark the centre of a compact cluster of massive stars at a distance of ∼12.5 kpc, believed to have hosted the magnetar progenitor star <cit.>. The physical size of the ring at this distance is ∼ 2.18×1.15 pc, it has a temperature ∼ 80 K <cit.>, and a flux of 0.4±0.1 Jy and 1.2±0.2 Jy at 16 and 24μm, respectively <cit.>. However, the inevitable difficulties in the analysis of this faint and complicated structure make these values rather uncertain and preclude a detailed assessment of the true ring morphology. The ring-like structure has been interpreted as due to illumination from nearby stars of a dust free cavity produced by the Giant Flare <cit.>. 


In this work we show the results of a series of radiation transfer (RT) calculations performed with a 3D dust radiative transfer code DART-Ray  <cit.> aimed at reproducing the emission of the infrared ring around SGR1900+14, assuming different plausible distributions for the dust illuminated by the nearby stars. In <ref> we introduce our initial assumptions, method, and the RT calculations. Results and Discussion follow in <ref> and <ref>.











§ RADIATION TRANSFER CALCULATIONS



In order to set up radiation transfer models appropriate to reproduce the observed ring emission, we first defined the distribution of stars and dust within the volume over which the RT calculations are performed (a cube of size 10 pc). In this section, we describe how this has been done given the constraints provided by the observations. Specifically, in <ref> the assumed ring distance, physical sizes and fluxes are given. In <ref> we explain how the dust distributions have been defined. In <ref> we describe how we found the appropriate viewing angles for the observer such that the projected ellipse is approximately of the same size and orientation as the dust emission ring. In <ref> we show how we derived the intrinsic 3D positions of the stars, relative to the magnetar, by using the constraint given by their projected positions on the sky. Finally, in <ref> we describe the specific RT models we considered and how the RT and dust emission calculations have been performed.    



 §.§ Ring distance, sizes and mid-infrared fluxes



We set up the size and geometry of the assumed 3D dust distributions by comparison with the observations (see Figure <ref>). We considered the distance measurement of <cit.>, which found d=12.5 ±1.4 kpc by using optical spectroscopy of the close-by stars. Assuming d=12.5 kpc, the lengths of the projected ring major and minor semi-axis are 2.18 and 1.15 pc. The major axis of the ring is rotated by about 22 from the R.A. axis. We considered the integrated fluxes for the dust ring in the mid-infrared as measured by <cit.>: F(24μm)=1.2±0.2 Jy and F(16μm)= 0.4±0.1 Jy.   






  



 §.§ Assumed ellipsoidal dust distributions


The observed 2D ring-like emission morphology observed on the Spitzer data is compatible with being the projection of a 3D ellipsoidal structure. We modelled it with a thin ellipsoidal shell, a uniform dust distribution around an ellipsoidal cavity or a more complicated distribution resulting from a stellar wind dust density profile which has been internally depleted of dust. In this section, we show how we have defined dust density profiles representative of each of these cases. The three profiles we have chosen should be considered as simplified representations of the complex dust distributions determined by the physical processes plausibly giving rise to the 2D ring we observe. Given the large uncertainties on the details of these physical processes, we could only assume simple shapes which qualitatively reproduce the expected dust distributions. 

The ellipsoidal surfaces, needed to define the above ellipsoidal structures, are described by the following formula for a constant normalized radius R:


    R^2=x^2/a^2+y^2/b^2+z^2/c^2


where a, b and c are the lengths of the three semi-axis of the "reference ellipsoid" with R=1. If for a given (x,y,z) we have R>1 or < 1, the point (x,y,z) belongs to an ellipsoidal surface which is either inside or outside the reference ellipsoid (but it has the same axis ratio). In order to consider the volume within an ellipsoidal shell, that is, the volume embedded between two ellipsoidal surfaces with different normalized radii, we consider all the points (x,y,z) satisfying the following relation:


    √(| R^2 -1| ) < Δ R
   

where Δ R represents the semi-width of the ellipsoidal shell (in normalized R units). For the elliptical shell distribution, the dust density is assumed to be constant within the volume defined by Eq.<ref> and zero outside. For the ellipsoidal cavity distribution, we assumed that the dust density is constant for R>1 and zero for R<1. Finally, for the stellar wind distribution, we assumed the following dust density radial profile: 

    ρ_d(R)=ρ_d(R_d)(R/R_d)^2     if R<R_d, 
    ρ_d(R)=ρ_d(R_d)(R/R_d)^-2     if R>R_d
 
with R_d=1. The wind density profile for R>1 decreases as R^-2 and thus resembles that expected in a stellar wind with elliptical symmetry. The profile for R<1 is hard to predict theoretically, since it is due to dust destruction processes with unknown parameters. For the sake of simplicity, we assumed it rises as R^2 until R=1. The three types of density profiles we defined are shown in Figure <ref>. 




 §.§ Derivation of the lines-of-sight reproducing the ring apparent sizes


For all the above dust distribution profiles,  the axis lengths and orientation of the ring are determined by the intrinsic lengths of the ellipsoid axis as well as by three angles: two angles (θ_ obs and ϕ_ obs) which specify the line-of-sight direction (see left panel of Fig. <ref>), and one angle that specifies the orientation of the 2D reference frame on the observer plane. By "observer plane", we mean the plane over which the ellipsoid is projected (that is, simply the plane of the data map). The ellipsoid projection is degenerate in the sense that different combinations of the ellipsoid intrinsic parameters can produce the same projected ellipse just by choosing appropriate observer viewing angles.

In order to be able to predict the shape of the projected ellipse for arbitrary combinations of the ellipsoid and observer parameters, we have used the formulae derived by <cit.>. Given an ellipsoid and the observer plane, these authors derived analytical formulae for the projected ellipse by determining all the points on the observer plane where the normal to the plane is tangent to the three-dimensional ellipsoid. By using equations 25 of GS81, we are able to predict the projected ellipse semi-axis and orientation given the ellipsoid semi-axis a, b, c and the observer line of sight angles θ_ obs and ϕ_ obs[The formulae in GS81 are written in terms of three rotation angles θ_GS, ϕ_GS and ψ_GS (see Fig. 3 in GS81) which can be connected to our definition of observer viewing angles θ_ obs and ϕ_ obs in the following way:

    θ_GS=π/2  ϕ_GS=π/2-θ_ obs  ψ_GS=ϕ_ obs-π/2.

By using the first relation, we force one axis of the 2D reference frame on the observer plane to be the projection of the z-axis of the 3D frame "xyz", where the ellipsoid is defined. In this way, θ_ obs and ϕ_obs are easily related to the angles ϕ_GS and ψ_GS in GS81 using the second and third relation above.]. 
Then, in order to handle the inverse problem of finding the combinations of the observer angles θ_ obs and ϕ_ obs that produce a projected ellipse with the same parameters of the observed dust emission ring, we wrote an optimization program which finds the right θ_ obs and ϕ_ obs for a given combination of ellipsoid semi-axis a, b, and c. Specifically, we choose to fix the values of b and c to 2.18 pc, the length of semi-major axis of the projected ellipse. We then assumed b/a= 2 or 4. For these different combinations of ellipsoid parameters, we found the values of θ_ obs and ϕ_ obs which allow to reproduce the shape and orientation of the projected ellipse (see <ref>). The values for θ_ obs and ϕ_ obs we derived are listed in Table <ref>. 





    


 §.§ Derivation of the intrinsic 3D positions of the stars from their sky position



In all the RT models that we calculated, we placed the magnetar at the origin of the 3D reference frame xyz, since the magnetar appears at the centre of the projected ellipse. For the two supergiant stars we had to find a way to place them such that their projected position on the observer plane coincides with their sky coordinates R.A. and Dec (see Table <ref>). 
Because this is the only constraint we have, each star can in principle be located at any point on a line perpendicular to the observer plane and intersecting the observer plane in (R.A., Dec), as shown in the right panel of Fig.<ref>. We parametrize the position of a star along this line in terms of its distance d_s from the magnetar. As shown in Fig.<ref>, for each value of d_s there can be up to two possible positions for the star (but also just one or even none if d_s is too small). If two possible positions exist, one of the two is chosen as location for the star. In Table <ref> we show the positions we derived for the stars assuming they are located at different distances d_s from the magnetar (that is, d_s = 1 and 3 pc). These distances are chosen such that the stars are either within the ellipsoidal cavity, or outside it but not too far from its border. We call these two types of star locations "IN" and "OUT" configurations.   




 











 §.§ 3D dust radiative transfer and dust emission calculations



In order to calculate the dust heating from the stars and the resulting dust emission, we used the 3D ray-tracing dust radiative transfer code DART-Ray   <cit.>. This code allows to solve the 3D dust radiative transfer problem for arbitrary distributions of dust and stars <cit.>: it follows the propagation of the light emitted by the stars within an RT model, including absorption and multiple scattering due to the dust. By calculating the variation of the light specific intensity I_λ throughout a RT model, it also derives the distribution of the radiation field energy density u_λ(x⃗) of the UV/optical/near-infrared radiation:

    u_λ(x⃗)=∫ I_λ(x⃗,n⃗)dΩ/c

 From u_λ and the dust density distribution, the dust emission can then be calculated at each position. We performed the dust emission calculation taking into account the stochastic heating of small grains <cit.>. This is important to consider since small grains tend to not achieve equilibrium with the interstellar radiation field and, heated by single photons carrying energies comparable to their internal energies, experience large temperature fluctuations. Thus, taking into account stochatical heating is important to correctly derive the dust emission spectra in particular in the mid-infrared. Unlike the equilibrium case, where a grain is characterized by only one dust temperature, in the stochastically heated case a grain has a certain probability P(T) to be at a certain temperature T. At each position inside the RT model, the probability function P(T) is calculated for each grain size a and composition k following <cit.> <cit.>. 
The function P(T) depends on both the local value of u_λ and the grain absorption coefficient Q^ abs_λ(a,k). Once P(T) is derived, the dust emission luminosity for a single grain ϵ_λ(a,k) can then be calculated as:

    ϵ_λ(a,k)=4π^2 a^2Q^ abs_λ(a,k)∫_0^∞ P(T)B_λ(T)dT
 
where B_λ is the Planck function. The total dust emission at each position is derived after integration over the grain size distribution (for each chemical species). 
Finally, by projecting the dust emission at each position onto the observer plane and by convolving the resulting maps to the instrument angular resolution (FWHM∼5.3 and 6 arcs for Spitzer 16 and 24 μm), dust emission maps can be derived as they would be obtained by an external observer. 

As input, DART-Ray  just needs a 3D Cartesian adaptive grid where for each element the dust density and the stellar volume emissivity are specified. Both distributed stellar emission and stellar point sources can be included. We created input grids representing the three types of dust distributions described in <ref>: an ellipsoidal shell, an ellipsoidal cavity and a ellipsoidal wind.
As we already mentioned in <ref>, we assumed  b=c=2.18 pc and b/a=2 or 4. In the case of the ellipsoidal shell dust distribution, we also assumed Δ R=0.10. We chose this value for Δ R because it produces a dust emission ring with approximately the same width of the observed ring. 

To set the dust optical properties (absorption and scattering efficiency, scattering phase angle parameter, dust-to-gas mass ratio), we assumed the dust model BARE-GR-S of <cit.> as implemented in the TRUST dust RT benchmark project <cit.>. In this dust model, the dust is composed by a mix of silicates, graphite and PAH molecules <cit.>. The abundance of each component in the dust model has been chosen to match observational data for the average extinction curves, chemical abundances and dust emission in the Milky Way.  

To set the dust/gas density in the RT model, we assumed several choices for the optical depth per unit parsec at 1 μm which, assuming the dust model mentioned above, correspond to a wide range of values for the hydrogen number density n_H. Since the density can vary within an RT model, the density values we will mention hereafter refer to the density at the reference ellipsoid (that is, for R = 1). Because of the long time required for a single RT/dust emission calculation (several hours), the density values have not been varied automatically in order to minimize the disagreement with the data but have been simply chosen in order to represent a typical diffuse interstellar medium (ISM) density (n_H=10 cm^-3) as well as much denser medium (n_H=100 and 1000 cm^-3, see Table <ref>). The ring integrated fluxes for the models with n_H=1000 cm^-3 are compatible with the observed fluxes, as it will be shown in the next section.

Following the procedure described in <ref>, we positioned the two M supergiant stars at either d_s=1 or d_s=3 pc from the magnetar (the two stars are at same distance from the magnetar but at two different 3D positions since they are projected in two different points on the observer plane). We assume these stars are emitting as blackbodies with effective temperatures and bolometric luminosities as found by <cit.> (see Table <ref>). 







§ RESULTS



In this section we show the results we obtained for the various assumed dust ellipsoidal distributions (see Figure  <ref>), dust densities (defined by the 1 μm optical depth per unit length or, equivalently, the hydrogen number density n_H), the stellar positions (which can be in the IN or OUT configurations) and the reference ellipsoid axis-ratio b/a. 
The parameters of the calculated RT models are shown in Table <ref>. In the same table we also show the fluxes we measured for the dust emission ring appearing on the 16 and 24 μm model maps. The ring photometry has been taken on an elliptical aperture containing the ring, in analogy to the measurement of <cit.> on the Spitzer maps. The 16 μm emission model maps, corresponding to all the RT models listed in Table <ref>, are shown in Figures  <ref>, <ref> and <ref>. The variation of the parameters listed above has a clear effect on the morphology of the synthetic maps and/or on the total flux. 










 §.§ Ellipsoidal shell models


The first two panels on the left in figure <ref> show the maps obtained for the ellipsoidal shell RT model with n_H=10cm^-3 and b/a=2 for the IN and OUT configurations.  As one can see, only the IN configuration gives rise to a clear dust ring shape, while the OUT configuration presents an additional bright feature within the ring–like shape. Although this enhancement resembles the one seen on the Spitzer maps, its nature is very different. The enhancement on the data maps appears at the position of the supergiant stars, it is a PSF-convolved point source and it is much brighter than the emission from the ring. Being much brighter than the expected blackbody emission for those stars, this emission is probably due to dust in the circumstellar material around the M supergiants. On the other hand, the bright feature observed in the OUT models is extended and shifted with respect to the position of the stars. This emission is due to dust in the elliptical shell that is closer to the stars located outside the shell and, thus, heated more strongly by them. 

The effect of modifying the b/a ratio, but keeping the stars inside the cavity, can be seen from the middle panel. Also in this case, although a dust ring is visible, an emission enhancement with comparable brightness to that of the ring is visible on the map. This is due to the more elongated shape of the cavity where some parts of the internal boundary of the cavity are  significantly more heated by the stars, thus causing this effect. In terms of emission flux, in this model the integrated ring flux is about a half of that of the model with b/a=2 (although the flux integrated over the entire map is higher). More importantly, the change of total flux due to the variation of b/a is too small to allow to reproduce the observed fluxes while keeping n_H=10cm^-3. This is because the b/a ratio can be varied, realistically, only within a rather small range of values (≈ 1-10).   

On the other hand, the dust density can in principle be varied over a wide range of values and, thus, can give rise to a substantial change in the total flux without affecting the emission morphology. We performed calculations with dust densities corresponding to n_H=100 and 1000 cm^-3 for a model with b/a=2 and stars within the cavity. As one can see from the right panels of Figure <ref> and the ring fluxes in Table <ref>, the predicted emission has the same morphology of the n_H=10 cm^-3 model but it is much more luminous. The dust emission SEDs of the ring for these two models are shown in Fig. <ref> (left panel). This plot shows that for the model with n_H=1000 cm^-3 the fluxes we obtained are within 1.5σ from those measured by <cit.>. 





 §.§ Ellipsoidal cavity models


We also performed RT calculations assuming a dust cavity geometry, where the dust is uniformly distributed outside the reference ellipsoid. This is expected in case the flare from the magnetar simply caused dust destruction within the dust cavity and did not affect significantly the density of the surrounding medium (thus, there is no significant enhancement of the gas/dust density close to the border of the cavity in this case). For n_H=10 cm^-3 and b/a =2, we calculated the maps for the IN and OUT geometry. The corresponding 16μm maps are shown in Figure <ref>. In both cases the emission is much more extended compared to that of the shell models (which do not contain any dust apart from that within the shell). However, the ring emission morphology is clearly recovered only in the case where the stars are inside the cavity. If the stars are outside, they illuminate the dust close to them very strongly. The emission from this dust dominates the emission seen on the maps. Instead, the ring–like emission is barely noticeable, with brightness close to that of the background emission. We point out that the bright dust emission seen in this case is much more extended than the point source emission seen on the data map at the position of the stars. We also calculated higher density models (n_H=100-1000 cm^-3) for the IN configuration. The maps shown in Fig.<ref> are characterized by a higher surface brightness, compared to the calculation with n_H=10 cm^-3, but a similar emission morphology. 

The ring fluxes for the ellipsoidal cavity model and for the IN configuration are listed in Table <ref>. These can be compared with the ring fluxes for the corresponding shell models with same parameters. As one can see, there are some differences between the shell and cavity models with same parameters, but the fluxes are similar. We also show the ring dust emission SEDs for the highest density models in Fig. <ref> (middle panel).  The cavity model with  n_H=1000 cm^-3 fits the observed integrated dust emission within the data error bars. This shows that both a thin dust shell and a dust cavity RT model give similar results for the ring integrated fluxes and morphology in the case that the stars are located inside the cavity. Also, overall, the results we obtained from the dust cavity models further suggest that the presence of the stars (or other radiation sources) within the cavity is a necessary requirement to recover the observed dust emission morphology. 



 §.§ Ellipsoidal dusty wind models


The last geometry we explored for the dust is that expected in the case the dust around the magnetar is distributed as in a stellar wind (with elliptical symmetry) and the dust has been only partially destroyed within the cavity region. The assumed density profile rises until the border of the cavity and then decreases as R^-2 as shown in Sect.<ref>. Given the previous results for the shell and cavity model, we have run only RT calculations assuming n_H=1000 cm^-3 at R=1 for this model. The results are shown in Fig. <ref> for both the IN and OUT configurations for the stars. In the case the stars are located inside, the morphology of the ring is recovered although not as clearly as we found before for some of the shell and cavity models.  There is substantial diffuse emission coming from the regions enclosed by the ring. Interestingly, a quite narrow and bright dust emission appears at the position of the stars. This feature resembles the bright dust emission source seen on the data at the star positions. A sharper profile, rising more rapidly close to R=1, would certainly reduce the brightness of this diffuse emission, as observed for the two previous configurations. However, given the large uncertainties on dust grain sizes or the properties of the destructing flare, we assume this dust density shape as a toy model to test the wind scenario.


We note that also for the wind model we do not recover the ring morphology in the case the stars are located outside the cavity region. The dust emission SED of the ring predicted in the case the stars are located inside the cavity region is shown in Fig.<ref> (right panel). The model is able to fit the total fluxes at 16 and 24μm (see also integrated fluxes in Table <ref>). This result confirms that a density of n_H∼1000  cm^-3 at the position of the border of the cavity region is necessary to fit the data, independently of the assumed dust distribution. 








 §.§ Comparison of the average surface brightness profiles for the high density models


In order to gain more quantitative insight on the similarities in morphology between the ring-like emission on the data and on the model maps, we compared average surface brightness profiles derived from the maps. We limited this comparison with the data to the models with the highest density (that is, n_H=1000 cm^-3), which are the only ones able to reproduce the total ring emission flux (see above). We derived the average surface brightness profiles in the following way. Firstly, we masked the emission from the brightest stars on the data maps and the corresponding regions on the model maps. Then, from the background-subtracted maps, we derived average profiles by averaging the values of the non-masked pixels within elliptical rings with same centre, axis-ratio and orientation as the infrared ring on the data. Finally, we normalized all profiles to their maximum value. The surface brightness profiles so derived at 16 and 24 μm are shown in Figure <ref> for the data and the high density RT models for all the assumed dust density distributions. The radial distance R_ map on the x-axis corresponds to the length of the semi-major axis of the elliptical rings used to derive the profiles. 

From both the 16 and 24 μm profiles several interesting features are evident. For radii larger than ∼1.2 pc, the profile for the shell dust distribution ends too sharply and it is thus unable to reproduce the tail of the emission observed on the data, which is much more extended. On the other hand, the tails of the profiles for the cavity and wind dust distributions decrease in a smoother way that is much closer to that observed on the data (see also average discrepancies in table <ref>). This finding implies that a thin dust shell is not able to reproduce the ring emission profile. Thus, the presence of a large amount of dust beyond the outer edge of the cavity is required. Furthermore, the comparison for the inner part of the emission profiles (R≤1.2 pc) shows that the wind dust distribution gives rise to excess emission inside the ring, which is incompatible with the data. A sharp rise of the dust distribution at small radii provides much better agreement with the data (see table <ref>). The preferred dust distribution is of an almost dust-free cavity with a sharp transition to a dust rich environment. Beyond the cavity the dust and/or interstellar medium densities decrease as a function of distance slower than the R^-2 wind-model.



  







§ DISCUSSION




We have performed several dust RT calculations assuming elliptical
dust shell/cavity geometries as well as a disrupted wind profile, and
by positioning the two supergiants stars inside or outside the dust
cavity. We have found that the dust ring morphology, similar to that
found on the Spitzer data of SGR1900+14, is recovered only in the
cases where the stars are inside the cavity. Furthermore, we
approximately reproduce the total integrated fluxes at 16 and
24 μm only by assuming a gas density of
n_H∼1000 cm^-3 for all the dust geometries we assumed.
The corresponding mass of the dust responsible for the ring emission
is M_ dust∼2 M_⊙. 


Given these results, the first question to ask is whether or not the
models that reproduce the observed dust emission morphology and total
flux are realistic or not. In particular the gas density, implied by
our modelling to explain the Spitzer infrared luminosity by dust
illumination, appears to be very high compared to that of the diffuse
galactic ISM.  Before discussing the possible nature of this high
density, we first clarify what assumptions/parameters in our modelling
might have caused an artificial high gas density, not representative
of the real ISM density around the magnetar. Firstly, we point out
that the gas density is not measured directly from the gas emission
but inferred from the dust density divided by a dust-to gas ratio of
0.00619 (which is characteristic of the assumed Milky Way dust model,
see section <ref>). However, this dust-to-gas ratio,
representative of the kpc scale ISM of the nearby Milky Way regions,
presents significant local variations in the ISM <cit.>. Furthermore, since the assumed size distribution
of the grains is also representative of the local Milky Way, this also
has an effect on the derived dust density. In fact, the MIR emission
in our modelling is mainly produced by small grains (sizes
a∼10^-3-10^-2 μm) which are stochastically heated. If the
grain size distribution is more skewed towards smaller grain sizes,
compared to the one we are assuming, this would require significantly
less dust mass to reproduce the observed MIR emission. The grain size
distribution is known to be affected by both dust destruction and
formation processes, but it is not possible to constrain it further
with our observations. On the other hand, we also note that if the
cavity has been created by dust destruction, the grain size
distribution there should instead favour the presence of large dust
grains rather than small ones <cit.>.  In
  fact, a number of studies <cit.> have
  shown that the X-ray flux is more effective at destroying small
  grains than larger ones. The precise evolution of the dust grain
  distribution is dependent on both the spectral shape and overall
  intensity of the illuminating source, on the composition of the
  grains, as well as on the relative importance of the processes of
  X-ray Heating, Coulomb Explosion and Ion Field Emission, the last
  two of which being particularly uncertain <cit.>. However, even within these uncertainties, all the models
  generally predict that smaller grains will be destroyed to larger
  distances than larger ones[If Coulomb Explosion dominates,
    then the distance to which grains of size a are destroyed by a
    source with energy spectral index α scales as
    a^-0.5-α/3 if grain charging is limited by internal
    energy losses, or as a^-1-α/2 if limited by the
    electrostatic potential.]. Therefore,  there would be a region in which only selective destruction took place, leaving behind a dust distribution skewed towards large grains. At the inner edge the distribution would be skewed towards big grains, progressively changing into the undisturbed (pre-burst) distribution at larger distances.     
  An attempt at modeling these
  effects would be worthwhile if the quality of the data were to allow
  a comparison with observations, but this is not possible with the
  current data.

Finally, in our modelling we only assumed the two supergiant stars,
the most luminous stars in the field, to be heating the dust
shell/cavity. However, other sources of radiation might well play a
role (e.g. other fainter stars within the cavity) and, in this case,
the needed gas density to match the observed fluxes would be lower. On
the other hand, note that the constant ∼10^34 erg/s X-ray
luminosity emitted by the magnetar is too low to power the dust
emission. In fact, the wavelength–integrated dust emission luminosity
for the models that fit the MIR fluxes is in the range 3.7-4 ×
10^35 erg/s. An additional mechanism to heat the dust is also
collisional heating in hot plasma, where the dust is heated by the
collisions with high energy electrons. This is expected if the dust is
embedded in shocked gas with temperatures of order of 10^6
K. However, in this case we might expect to see an X-ray diffuse
emission from the hot gas around the magnetar as well, which is not
observed.



Vrba et al.(2000) argued that SGR 1900+14 is associated with a cluster
of young stars (much fainter in apparent magnitude than the two M
supergiants we considered in our modelling) which are probably
embedded in a dense medium. This interpretation is qualitatively
consistent with our results.  As proposed by <cit.>, the
1998 Giant Flare could have produced the cavity by destroying the dust
within it. Assuming a constant dust density within the cavity region,
corresponding to n_H=1000 cm^-3, we estimated a total dust mass
of order of 3 M_⊙ that was plausibly present before being
destroyed by the flare. An energy of about E∼ 6×10^45erg
would suffice to destroy this amount of dust, consistent with the
estimates by <cit.> based on Eq. 25 in
<cit.>.  The size of the region with destroyed dust
would be larger for smaller grains, as discussed above.
In this scenario, the high density we
derived would be similar to the high density ISM around the
magnetar. Furthermore, we note that high density of the ISM
(n_H=10^5–10^7  cm^-3) has been found in the environment
surrounding GRBs <cit.>, which should be similar to that
where magnetars are located.

The wind model we considered was meant to be similar to the scenario where the dust distribution outside the cavity was mainly determined by the wind of the magnetar progenitor while internally disrupted by the Giant Flare. However, gas densities at 1pc distance in a typical stellar wind are expected to be several orders of magnitudes lower than those we found (n_H∼1 cm^-3,  estimated from the mass-loss rates Ṁ of <cit.> assuming n_H ∝Ṁ/v/r^2 with v = 1000 km/s). Thus, this last scenario is unlikely if the ring density is indeed so high.  


 
Another possibility is that the dust emission ring is the infrared emission from the supernova remnant (SNR) of the magnetar progenitor. The dust mass associated with the shell model with n_H=1000 cm^-3 is M_dust=1.9 M_⊙, which is a factor 3–4 higher than the measurement of dust mass around SN1987 by <cit.>. However, given the large uncertainties in the inferred dust masses, and that our value for the dust density might have been overestimated because of the reasons given above, it may well be that the amount of dust needed for the shell model is compatible with that of SNRs. The SNR scenario was also considered by <cit.> but discarded because of the lack of observed radio and X-ray emission from the ring. However, if we consider i) the IR/X luminosity ratio of ∼10^-1-10^2 measured by <cit.> for many SNRs, ii) the total IR luminosity of the ring in our models (∼4×10^35 erg/s), and the X-ray detection limit for the ring <cit.>, this structure is still compatible with being a SNR with high IR/X ratio.

On the other hand, we can also compare the 24μm luminosity with the expected X-ray luminosity, according to Figure 12 of <cit.>, which studied a sample of SNR in the Large Magellanic Cloud. If the magnetar was located at the distance of the LMC (50 kpc), we would have ν F_ν(24μ m)= 1.5×10^-10 erg/s/cm^2, and the relative expected X-ray flux would be 2×10^-10 erg/s/cm^2. This would translate in an intrinsic X-ray luminosity of ∼6×10^37 erg/s, which should have been clearly detected in the case of the SGR 1900+14. Hence, given the information we have at hand we cannot discard the SNR scenario on the base of the observed IR/X-ray luminosity ratio, although it would be a rather peculiar remnant compared with what we see around other Galactic pulsars or magnetars <cit.>. 

However, if we also consider the shape of the normalized average surface brightness profiles, shown in Fig.<ref>, this provides a strong evidence that 1) there is very little amount of dust inside the cavity and 2) the emitting dust is much more extended than a simple thin shell. These findings are compatible with the scenario where the cavity has been produced by the Giant Flare within a high density medium. However, the SNR scenario would still be acceptable in the case the transition in density between the shell and the surrounding ISM is smoother than what we assumed in our modelling. 

Regardless of the origin or the exact distribution of the illuminated dust, or the exact nature of the dust free cavity, our models show that we are able to observe this illuminated dust structure only because of two favourable characteristics: 1) the high dust density in the local region, and 2) the illuminating stars coincidentally lay inside the shell. Similar dust structures might potentially be present around many other magnetars or pulsars but they would be invisible to us because of the lack of either one of the two above local properties of this particular object.

 


GN acknowledges support by the EU COST Action MP1304, for the Short Term Scientific Mission where this project was completed, and by the Leverhulme Trust research project grant RPG-2013-41. N.R. acknowledges funding in the framework of the NWO Vidi award A.2320.0076, and via the European COST Action
MP1304 (NewCOMPSTAR). N.R. and D.F.T. are supported by grants AYA2015-71042-P and SGR2014-1073.  RP acknowledges support from the NSF under grant AST-1616157. JMG is supported by grant AYA2014-57369-C3-1-P.




























[Camps et al. (2015)]Camps15
Camps, P. et al., 2015, A&A, 580, 87  

[Davies et al. (2009)]Davies09
Davies B. et al. 2009, ApJ, 707, 844  

[Draine (2003)]Draine03
Draine, B. T. 2003, ARA&A, 41, 241  

[Duncan & Thompson (1992)]Duncan92
Duncan, R. C. & Thompson, C. 1992, ApJ, 392, 9  

[Eikenberry et al. (2001)]Eikenberry01
Eikenberry, S. S. et al. 2001, ApJ, 563, 133  

[Fruchter et al. (2001)]Fruchter01
Fruchter, A. et al. 2001, ApJ, 563, 597  

[Gaensler et al. (2001)]Gaensler01
Gaensler, B. M. et al. 2001, ApJ, 559, 963  

[Gendzwill & Stauffer (1981)]Gendzwill81
Gendzwill D.J. & Stauffer M. R. 1981, Journal of the International Association for Mathematical Geology, 
Volume 13, Issue 2, pp 135-152  


[Green (1984)]Green84 
Green D. A., 1984, MNRAS, 209, 449 


[Hurley et al. (1999)]Hurley99 
Hurley K., et al., 1999, Nature, 397, 41 

[Hurley et al. (2005)]Hurley05
Hurley, K. et al. 2005, Nature, 434, 1098  

[Israel et al. (2008)]Israel08
Israel, G. L. et al. 2008, ApJ, 685, 1114  

[Kaplan et al. (2003)]Kaplan03
Kaplan D. L. et al., 2003, ApJ, 590, 1008  

[Koo et al. (2016)]Koo16
Koo, B. et al. 2016, ApJ, 821, 20   

[Kudritzki (2002)]Kudritzki02
Kudritzki, Rolf P. 2002, ApJ, 577, 389  

[Lazzati & Perna (2002)]Lazzati02
Lazzati, D., Perna, R.,  2002, MNRAS.330, 383

[Martin et al. (2014)]Martin14 
Martin J., Rea N., Torres D.  F., Papitto A., 2014, MNRAS, 444, 2910 

[Matsuura et al. (2011)]Matsuura11
Matsuura, M. et al. 2011, Sci, 333, 1258    

[Mazets et al. (1979)]Mazets79 
Mazets E. P., Golentskii S. V., Ilinskii V. N., Aptekar R. L., Guryan I. A., 1979, Natur, 282, 587 


[Muno et al. (2006)]Muno06 
Muno M. P., Law C., Clark J. S., Dougherty S. M., de Grijs R., Portegies Zwart S., Yusef-Zadeh F., 2006, ApJ, 650, 203 


[Natale et al. (2014)]Natale14
Natale G. et al., 2014, MNRAS, 438, 3137 

[Natale et al. (2015)]Natale15
Natale G. et al., 2015, MNRAS, 449, 243 

[Olausen & Kaspi (2014)]Olausen14
Olausen, S. A. & Kaspi, V. M. 2014, ApJS, 212, 60    


[Perna & Lazzati (2002)]Perna02
Perna, R. & Lazzati, D. 2002, ApJ, 580, 261  

[Perna, Lazzati & Fiore (2003)]Perna03
Perna, R., Lazzati, D. & Fiore, F. 2003, ApJ, 585, 775 

[Perna & Pons (2011)]Perna11
Perna, R. & Pons, J. A., 2011, ApJ, 727L, 51

[Rea & Esposito (2011)]Rea11
Rea, N. & Esposito, P. 2011, ASSP, 21, 247  

[Reach et al. (2015)]Reach15
Reach, W. T. et al. 2015, ApJ, 811, 118

[Seok et al. (2013)]Seok13
Seok, J. Y. et al. 2013, ApJ, 779, 134    

[Steinacker et al. (2013)]Steinacker13
Steinacker, J. et al. 2013, ARA&A, 51, 63


[Thompson & Duncan (1993)]Thompson93
Thompson, C. & Duncan, R. C. 1993, ApJ, 408, 194 

[Thompson & Duncan (1995)]Thompson95
Thompson, C. & Duncan, R. C. 1995, MNRAS, 175, 255  

[Turolla et al. (2015)]Turolla15
Turolla, R. et al. 2015, RPPh, 78, 6901   


[Voit (1991)]Voit91
Voit, G. M. 1991, ApJ, 379, 122  

[Vrba et al. (1996)]Vrba96
Vrba, F. J. et al. 1996, ApJ, 468, 225  

[Vrba et al. (2000)]Vrba00
Vrba, F. J. et al. 2000, ApJ, 533, 17   

[Wachter et al. (2008)]Wachter08
Wachter et al. 2008, Nature, 453, 626   

[Waxman & Draine (2000)]Waxman00
Waxman, E. &  Draine, B. T. 2000, ApJ, 537, 796   


[Zubko et al. (2004)]Zubko04
Zubko et al., 2004, ApJS, 152, 211   












