Mon.Not.R.Astron.Soc.000,1–16(2004) Printed2February2008 (MNLATEXstylefilev2.2) Simulations of spectral lines from an eccentric precessing accretion disc Stephen B. Foulkes1, Carole A. Haswell1, James R. Murray2,3, Daniel J. Rolfe 1,2 1Department of Physics & Astronomy, The Open University,Walton Hall, Milton Keynes,MK7 6AA, UK. 4 2Department of Physics & Astronomy, Universityof Leicester, UniversityRoad, Leicester,LE1 7RH, UK. 0 3Department of Astrophysics & Supercomputing, Swinbourne Universityof Technology, Hawthorn, VIC 3122, Australia. 0 Email SBF: [email protected], CAH: [email protected], JRM: [email protected], DJR: [email protected] 2 n a Accepted. Received J 9 1 ABSTRACT Two dimensional SPH simulations of a precessing accretion disc in a q = 0.1 binary 1 system (such as XTE J1118+480) reveal complex and continuously varying shape, v kinematics, and dissipation. The stream-disc impact region and disc spiral density 0 wavesareprominentsourcesofenergydissipation.Thedissipatedenergyismodulated 7 ontheperiodP =(P −1−P −1)−1withwhichtheorientationofthediscrelative 3 sh orb prec to the mass donor repeats. This superhump modulation in dissipation energy has a 1 variation in amplitude of ∼10% relative to the total dissipation energy and evolves, 0 4 repeatingexactlyonlyafterafulldiscprecessioncycle.Asharpcomponentinthelight 0 curveisassociatedwithcentrifugallyexpelledmaterialfallingbackandimpactingthe / disc. Synthetic trailed spectrograms reveal two distinct ‘S-wave’ features, produced h respectivelybythestreamgasandthediscgasatthestream-discimpactshock.These p S-waves are non-sinusoidal, and evolve with disc precession phase. We identify the - o spiral density wave emission in the trailed spectrogram. Instantaneous Doppler maps r showhowthestreamimpactmovesinvelocityspaceduringanorbit.Inourmaximum t s entropy Doppler tomogram the stream impact region emission is distorted, and the a spiral density wave emission is suppressed. A significant radial velocity modulation : v of the whole line profile occurs on the disc precession period. We compare our SPH i simulation with a simple 3D model: the former is appropriate for comparison with X emissionlineswhilethelatterispreferableforskewedabsorptionlinesfromprecessing r discs. a Key words: accretion: accretion discs - X-rays: binaries - binaries: low mass X-ray binaries, cataclysmic variables - stars: individual: A0620-00, XTE J1118+480,GRO J1655-40- methods: numerical 1 INTRODUCTION modulated tidal stresses (Osaki1985; Murray 1998). This modulation is accompanied by a change in the solid angle Superhumps are photometric modulations in the observed subtended by the disc at the central X-ray source, and it lightcurvesfoundinsomeshortperiodcataclysmicvariables is this shape change which leads to the irradiation-powered (CVs) (Patterson 1998), and in some of the transient low superhumpsin X-ray bright SXTs (Haswell et al. 2001). mass X-ray binaries known as soft X-ray transients (SXTs) MaterialfromthedonorstarleavesthefirstLagrangian (Zurita et al. 2002). In general the superhump period, P , sh point(L1)andflowsalongaballistictrajectorytowardsthe isafewpercentlongerthanthebinaryorbitalperiod,P , orb centralprimaryobject.Whenthisgasstreamencountersthe andisattributedtotheprogradeapsidalprecessionofanon- accretion disc material, large amounts of energy are dissi- axisymmetric accretion disc. Such precession is expected in patedinabrightspotregionontheouteredgeofthedisc.As extreme mass ratio systems (q . 0.3, where q = M2/M1 the eccentric accretion disc precesses in the inertial frame, is the mass ratio) in which the disc encompasses the 3:1 and the orbital motion of the binary proceeds, the energy resonance (Lubow 1991). If the disc extends to the 3:1 res- dissipated at thebrightspot will bemodulated as described onant radius, it becomes eccentric and non-axisymmetric byRolfe, Haswell & Patterson (2001 hereafter RHP2001). (Murray 1998)anditwillthenprecessunderthetidalinflu- enceofthesecondarystar.InCVs,superhumpsinthelight Spectroscopic observations of superhumping CVs curve are understood as being powered directly from the (Vogt 1982; Hessman et al. 1992; Patterson et al. 1993a) 2 S. B. Foulkes et al. and SXTs (in particular XTEJ1118+480) have shown that The model includes a simple 3D structure for the trailedspectrogramsofsuchobjectsarecomplexandchange brightspotregion.Anellipticalcross-sectioninther-zplane from night to night (Torres et al. 2002; Zurita et al. 2002; is assumed, with size r and h 1 centred on the outer spot spot Haswell et al. 2004).Overintervalsofdaysormonths(many discedge.r ,h andthespotbrightnessdecreasedown- spot spot times the orbital period) the spectral line profiles evolve streamfrom theinitialimpactase−θ2/∆θs2pot whereθ isthe from almost symmetric to enhanced red-shifted or blue- angular distance downstream from the initial impact point. shifted. These changes have been attributed to a precess- Upstream of the impact the brightspot surface is rounded ingaccretion disc,andhereinwepresentcalculations which offwithahemisphereofuniformbrightnessequaltothatat strongly support this interpretation. theinitialimpact.Weusedr =0.03aandh =0.012a spot spot In this paper we describe the results from two simula- at θ=0 (aisthebinaryseparation). Theangularextentof tions of a precessing eccentric accretion disc. The first sim- thebrightspotregion isset by∆θspot =∆θ0adisc/r0,where ulation consisted of a simple analytical three-dimensional r0 is the disc radius at the impact point and adisc is the representation of the disc and brightspot region based semi-major axis of the disc. This keeps the arc length of on RHP2001. The second model was a high-resolution the brightspot region roughly constant at ∆θ0adisc as the smoothed particle hydrodynamics (SPH) simulation. SPH eccentric disc precesses. Weused ∆θ0=30◦. is a Lagrangian numerical method for modelling fluid flow The disc was modelled as an ellipse with a =0.31a disc using a set of moving particles. The fluid properties at any and eccentricity 0.093 at the outer edge2. Streamlines and point within the flow field are determined by interpolating brightnesscontourswithinthediscareellipseswithsizeand from the local particle properties. This interpolation takes eccentricity decreasing smoothly tozero at theprimary ob- the form of a summation over the local particles, weighted ject,withthevelocitygivenby|V~ |.Theinnerboundaryis disc accordingtotheirdistancefrom theevaluation point.Fora near-circular with radius 0.025a. The disc brightness varies comprehensive review of the method see Monaghan (1992). as a−1.2. The disc is flared with total thickness increasing orbit In section 2 we detail the binary system modelled and linearly from 0 at the primary object to 0.02a at the outer the two numerical methods employed to model the system. edge.Anorbitalinclinationof70◦wasusedtoavoideclipses. Section 3 presents our results and Section 4 is devoted to In the brightspot region, where the stream and disc discussion and conclusions. flows merge, we adopted a velocity distribution which was non-zerobetweenthelocaldiscvelocityatthatpoint,V~ , disc and the velocity V~disc + (V~stream,0 − V~disc)e−θ2/∆θs2pot,v. 2 SIMULATIONS Thus, at the impact point, velocities between V~ and the disc 2.1 Binary Parameters ballistic stream velocity at theimpact, V~stream,0 contribute to the Doppler broadening of the local line profile. Down- TosimulateatypicalSXT,weadoptedamassratioq=0.1. stream the local velocity smoothly approaches the disc ve- Theorbitalperiod,P ,is 4π2a3/GM 1/2,whereaisthe locity. Upstream of the impact the exponential factor was orb t orbital separation, and Mt(cid:0)is the total (cid:1)mass of the binary set to 1. ∆θspot,v behaves as ∆θspot as the eccentric disc system. The mass transfer rate was 1.46×10−9Mtyr−1. precesses. Weused ∆θ0,v =5◦. The local linewidths were assumed to have Gaussian profiles, with the disc width set to the disc thermal veloc- 2.2 Simple Analytical Approximation to a ity. The brightspot profile width was the quadratic sum of Precessing Eccentric Disc thediscthermalvelocityandthethermalvelocityofa104K Herewebrieflydescribeasimplemodelfortheeccentricdisc sourcedecayingase−θ2/4∆θs2pot downstream.Thedecaywas and stream-disc impact developed from that used to model imposed for consistency with the brightspot brightness, as- latesuperhumpsintheSUUMaCVIYUma,seeRHP2001. suming brightspot brightness ∝Tb4rightspot. The model assumes the brightspot emission is directly pro- portional tothekineticenergy oftherelativemotion of the 2.3 SPH Simulation stream and disc flows at the impact point, i.e. the rate of energy dissipation at thebrightspot is: The accretion disc was also modelled using a two- dimensional SPH computer code that has been described E˙ = 21M˙Vr2el⊥+ 12M˙Vr2elk VVrelk in detail in Murray (1996, 1998). Modifications have been relk made to convert the code from FORTRAN to C++. C++ (cid:12) (cid:12) where M˙ is themass flowrate of thestr(cid:12)eam(cid:12)and V and is an Object Orientated computer language allowing each rel⊥ particle to have its own private attributes (e.g. position, V aretherelativevelocityofthestreamanddiscresolved relk velocity) and particle interactions can be coded in a more perpendicularandparalleltothediscrespectively.Thedisc flexiblemanner.Thecodewasalsoextendedtoallowitsex- velocity at radius r from the primary object is assumed to ecutiononanarrayofprocessorsusingtheMessage-Passing bethatofanellipticalorbitaroundtheprimaryobject,that Interface Standard3. The SPH results presented here were is: generated on a small private Linuxnetwork. |V~disc|=rGM1 2r − a 1 , In the SPH simulation the total system mass, Mt, and (cid:16) orbit(cid:17) where aorbit is the semi-major axis of the orbit and M1 is 1 SeeFigure8ofRHP2001. themassoftheprimarilyobject.Thevelocityofthestream 2 TheeccentricitywaschosetomatchthatofourSPHsimulation V~stream is simply that of theballistic trajectory. 3 http://www-unix.mcs.anl.gov/mpi/ Spectroscopic simulations of superhumps 3 Figure 1. The top row (a, b, c) are accretion disc surface density maps using a logarithmic colourscale for binary phases 0.00, 0.36, and0.64respectively. Thesecondary orbitsanti-clockwisewithrespecttothe inertialframe,withmassbeingadded fromthe L1 point attherightofeachmap.ThesolidcurveistheRochelobeoftheprimary,plottedinaframethatco-rotateswiththebinary.Thearrow in the lower right-hand corner of each map indicates a fixed direction to an observer. The middle row (d, e, f) are the corresponding dissipationmapsplottedwithalogarithmiccolourscale.Thebottomrow(g,h,i)arethematchingDopplermapsplottedusingalinear colourscale. The letter “G” on map i, indicates the point at which the gas stream impacts on the outer edge of the accretion disc and theletter “H”indicatesthearcofemissionproducedbythediscparticlesintheimpactregion. 4 S. B. Foulkes et al. the binary separation, a, were both scaled to unity. The occur and hence the hot spot region will have a reduced binaryorbital period, P ,was scaled to2π.Theaccretion contributiontotheenergydissipated.However,thehotspot orb disc had an open inner boundary condition in the form of region will still contribute a significant amount to thetotal a hole of radius r1 = 0.025a centred on the position of the dissipation observed from the accretion disc. primary object. Particles entering the hole were removed fromthesimulation.Particlesthatre-enteredthesecondary Roche lobe were also removed from the simulation as were 3 RESULTS particlesthatwereejectedfrom thediscandhadadistance >0.9a from thecentre of mass of theprimary. 3.1 Disc morphology: orbital modulation Theartificial viscosityterm ofMurray(1996) was used In this subsection we follow the evolution of the SPH sim- toimprovetheeffectiveshearviscosityintheSPHequations ulated accretion disc through a complete binary orbit. The of motion. In the inner regions (<∼0.12a) of the accretion images(a,b,c)inFigure1areasetofsurfacedensitymaps, disctheparticlesgenerallyexecutedcircularorbits,andthe plottedintheco-rotatingframe,fororbitalphases0.00,0.36 radialsheartermdominatedtheinnerdiscdissipation.How- and 0.64 respectively. Each map uses the same logarithmic ever,intheouterdiscthetidalinteraction ofthesecondary colourscale,withbluerepresentinglowdensityandredhigh, was stronger and the bulk term was more significant. The we have also plotted the primary Roche lobe. Since we can discenergydissipationwascalculatedusingtheSPHenergy view the disc from any angle, we have arbitrarily fixed the equationasinMurray(1996).TheShakura-Sunyaev(1973) conjunctionoftheprimary(andhencezeroorbitalphase)to viscosityparametersweresettoα of0.1andα of1.0, low high correspondtoFigure1(a).Thearrowinthelowerright-hand and theviscosity switched from low to high as described in cornerofeachmapindicatesafixeddirectionintheinertial Truss et al. (2000). The SPH smoothing length, h, was al- frame and we calculate radial velocities as they would be lowed to vary in both space and time and had a maximum observed by an observer in thisdirection. valueof 0.01a. These density maps show the disc to be asymmetric, The SPH code was started with no particles, i.e. zero eccentric and continuously changing shape under the tidal mass in the accretion disc. A single particle was injected into the simulation every 0.01Ω−1 at the first Lagrangian influence of the donor star; its overall position is moving orb only very slowly in the inertial frame. There are wrapped point. The sound speed in the donor atmosphere was used spiraldensitywavesthatextendfromtheoutermostregions astheinitialspeedofeach particleandtheyfollowed abal- ofthediscdowntoapproximatelythecircularisation radius listic trajectory. Within twenty binary orbits the resulting of the system. These produce shear and dissipation in the discbecameeccentricandencounteredtheLindblad3:1res- outer regions of the disc. The spiral density waves propa- onance.Afterthedischadreachedamassequilibriumstate, gate angular momentum outwards, thus allowing the disc wefollowed itsevolutionthroughfourdiscprecessioncycles gastomoveinwardtowardsthecentralprimaryobject.The (∼180 binary orbits) to produce the results described here. precession of the disc means that the density maps do not There were approximately fifty-three thousand particles in repeat from one orbit to another; instead the maps repeat thedisc,andtheaveragenumberof‘neighbours’,thatisthe on a period a few percent longer than the orbital period. average number of particles used in the SPH update equa- We have placed a Microsoft AVI movie of the change tions, was 8.9 particles. in surface density as function of orbital phase on an In- Thecomplicated hydrodynamicsofthestream-discim- ternet web-site (see http://physics.open.ac.uk/FHMR/ for pacthavebeeninvestigatedanalyticallybyHessman(1999) more information). and numerically by Armitage & Livio, first with the SPH technique (1996) and later with the grid based ZEUS code (1998). In the latter paper, they found the character of the 3.2 Light curve stream-overflow in particular depended strongly upon the equationofstate.Kunze,Speith&Hessman(2001)compre- InFigure1(d,e,f),wehaveplottedviscousdissipationmaps hensively studied the impact with much higher resolution correspondingtothedensitymapsinFigure1(a,b,c).Each SPH simulations, for a range of binary mass ratios. They dissipation map is represented using the same logarithmic showedthatasubstantialfraction(perhapsmorethanhalf) colourscalewithbluelowdissipationandredhigh.Itisclear of the stream overflowed the outer disc edge, to impact the from this set of images that the inner disc is in a state of disc close to the circularisation radius. The latter paper is constanthighdissipation.Thespiralcompressionwavesand particularly relevent as it demonstrates that the detailed the impact of the gas stream on the outer edge of the disc shock structure of the impact region can be captured us- also generate high levels of dissipation. ingtheSPHtechnique,ifsufficientnumbersofparticlesare We are unable to directly resolve the morphology of deployed. The computational cost is however so great as to the accretion disc in observed interacting binary star sys- make a calculation following the evolution of a disc pro- tems. The simplest directly observable quantity which re- hibitively expensive. veals some of the underlying properties is the light curve. As our SPH calculations are two dimensional, stream- Hence we use the SPH model to simulate the disc lumi- overflowcannotoccurandthebrightnessofthehotspotwill nosity and generate an artificial light curve. We assume beoveremphasisedbyafactordeterminedbycooling inthe thattheluminosityisgeneratedthroughviscousdissipation stream. The SPH results presented in the following section which is promptly radiated away from the point at which will show that the complex stream-disc interaction region it was generated. The simulated light curve can be directly dominates the energy dissipation signature of the accretion compared with the observations. We must, however, bear disc. In a three dimensional SPH code stream-overflow will in mind that in observed systems there will generally be Spectroscopic simulations of superhumps 5 an orbital modulation in the visibility of each location in theco-rotating frame. Hence, oursimulated light curves,in which we have essentially assumed 100% visibility for each location throughout theorbit, are likely to be simpler than observed light curves. The upper left plot of Figure 2 shows five simulated light curves. The upper curve corresponds to the total vis- cous dissipation from the whole disc; the lower four curves correspond to the four smaller annuli defined in Figure 2 caption.Thedissipationintheinnerdisc(.0.12a,interme- diatebetween‘c’and‘d’inFigure2)isconstant,exceptfor stochastic noise. The lower left-hand plot is the dissipation from the innermost disc region and was generated by sub- tractingdissipationlightcurvebfromdissipationlightcurve a.Subtractingcurvecfromaandsubtractingcurvecfromb would also produce similar constant dissipation curves. Su- perimposed on this, for disc regions (r > 0.2a) there is a repeating series of “humps.” These humps recur on the su- perhumpperiod,P .Themodulation consistsoftwosepa- sh ratecomponents.Thefirstistherelativelysmoothcontinu- ous modulation, which has a minimum at orbital phase 0.0 Figure 3. Power spectrum of the SPH dissipation light curve in Figure 2 (corresponding to Figure 1(a,d)) and reaches a Figure2(a).Thetimeseriescoveredapproximately60superhump maximum at orbital phase ∼0.64, (corresponding to Figure periodswitharesolutionof0.01Ω−or1b.Thepeaklabelled0isthe 1(c,f)). The second is the much sharper spiky signal with superhumpfrequency,thesuperhumpharmonicsarealsolabelled a short duration which appears at orbital phase ∼0.36 and (1-5).ThelowerpeakscorrespondtolinearcombinationsofΩorb, ω andΩ. corresponds to Figure 1(b,e). Photometric superhumps are invariably detected in optical light. Since the hot inner re- gionsofaviscously-heateddisccontributerelativelylittleto placedMicrosoftAVImoviesshowingthesurfacedissipation the total optical light, the comparison between observation in the SPH simulation and in the analytical model on the and our simulation might best be made by comparing one Internetsite http://physics.open.ac.uk/FHMR/ . of the lower curves in Figure 2. To determine the superhump period, P , we took a sh Figure 1(d) corresponds to the light curve minimum. Fourier transform of approximately sixty superhump cycles Fromthismapitisclearthatthedissipationatthestream- ofthesimulatedlightcurve.Thepowerspectrumisshownin discimpactisataminimum.Thisisbecausethespiralcom- Figure 3. The spectrum is plotted using a logarithmic scale pression arm in the outer disc has moved such that gas is torevealthelowerlevelpeaks.Alargepeakthatcorresponds escaping from theprimary Rochelobe and returningto the tothesuperhumpperiod,labelled0inthefigure,dominates secondary.ThegasleavingtheL1 pointfromthesecondary the spectrum; the superhump harmonics have also been la- encounters the disc gas instantly and the two flows have belled(1-5).Usingthepowerspectrumweestimatedthesu- a low relative velocity. Hence relatively little dissipation is perhumpperiodtobe(1+0.0222±0.0005)P ;thisimplies orb generated by the converging flows. At this phase the dissi- P =(46±1)P .Thesuperhumpperiodmanifestinour prec orb pation is dominated by the shear flow in the inner regions simulation is in agreement with the analysis of Mineshige, of the disc. Hirose&Osaki(1992).Theothermuchloweramplitudepe- Asthebinaryorbitproceeds,alargevoidopensupbe- riodicitiespresentoccuratfrequenciesΩrelatedtoΩ and orb tweentheL1pointandtheaccretiondisc,asshowninFigure disc apsidal precession frequency ω by: 1 (b) and (e) and the spiral arm which previously extended towardstheL1regionbecomesdetachedfromthediscunder k(Ω−ω)=j(Ω−Ωorb) centrifugalforcesandthenfallsbacktowardsthedisc.When where k,j are positive integers. thisgashitsthedisc,theimpactcausesasuddenriseindis- sipation, corresponding to thespikypeak in thelight curve atorbitalphase∼0.36,andtotheupperinhomogeneousarc 3.3 Trailed Spectra on Figure 1(e). Once this gas has been re-assimilated into thediscflowthedissipation falls again. Thesecondary con- Kinematics in observed systems are revealed primarily by tinuesitsorbitaroundthedisc,andthestreamflowandthe thespectral line profiles. To compare the results of the two disc flow at the stream-disc impact point become increas- simulations with CV and SXT observations we therefore inglymisaligned.Astherelativevelocityofthetwoflowsin- generated syntheticemission line profiles and used these to creases,thedissipation inthebrightspotregioncorrespond- build up synthetictrailed spectrograms. inglyincreases,untilitreachesamaximumatorbitalphase An accretion disc emission line is Doppler broadened. ∼0.64,showninFigure1(c)and(f).Astheorbitproceeds, Horne&Marsh(1986) calculated spectrallineprofilesfrom the angle between the stream and disc flows at the impact discswithsmoothazimuthallysymmetriclineemissivitydis- point decreases again, until the next minimum is reached tributions, and Horne (1995) further included the effects of at orbitalphase1.022. Similarbehaviourisseen inthesim- anisotropic turbulence on the local line broadening profile. ple3Dmodel,andwasfullyexploredinRHP2001.Wehave We generated artificial emission line spectra for emissivity 6 S. B. Foulkes et al. Figure2. ThetopleftplotisSPHviscousdissipationlightcurvesfordifferentregionsoftheaccretiondisc.Thetopcurve,labelleda,is forthewholeaccretiondisc,theneachcurveindescendingordercorrespondstodissipationfromthediscatradii>0.05a(b),>0.1a(c), >0.2a(d) and>0.3a(e)respectively. Thelightcurveminimumandtwomaximumpointsareindicated. Thelowerleft-handplotisthe signalfromthediscinnerregionandwasgeneratedbysubtractinglightcurvebfromlightcurvea.Theplotsontheright-handsideare discdissipationmapsthatcorrespondtothelightcurveminimumandtwomaxima.Thecircleslabelledb,c,dandeshowthedistances fromtheprimarilythatcorrespondtothelightcurvesb,c,d,erespectively. distributions calculated from both the simple 3D code and is shown in Figure 5 and the trailed spectra for a full disc the2D SPH code. precession cycle is shown in Figure 6. For the simple 3D model we adopted Gaussian local UsingtheSPHmodelthedissipationfromtheaccretion line profiles with widths proportional to the disc thermal discisdividedintoradialvelocitybins.Theemissionineach velocity; the brightspot emission was assumed to be ther- velocity bin is the sum of the viscous dissipation from each mally broadened using a brightspot temperature of 104K. particle within the velocity bin. For each SPH time step, We did not check for occultation of the accretion disc by 0.01Ω−1, an emission line profile was generated and shown orb the brightspot, but the 70◦ inclination was chosen to avoid in the trailed spectrogram in Figure 7. We have used a lin- eclipses by the secondary. Figure 4 shows three dissipation eargreyscale to representdissipation intensity.Thebottom mapsfromtheanalyticalmodelfororbitalphases0.00,0.36 of the spectrogram corresponds to the minimum flux point and 0.64 which correspond to the same orbital phases used on the artificial light curve of Figure 2 (labelled ‘min’). At in Figure 1. These maps show the dissipation plotted using approximately orbital phase 0.36 (‘max1’) there is a large alogarithmicgreyscalewithwhitelowdissipationandblack emission associated with the first large spiky hump on the high.TheprimaryRochelobehasbeenplottedusingasolid light curve. At about orbital phase 0.64 (‘max2’) there are curve.Thearrows indicatea fixeddirection, as in Figure 1. two distinctive emission arcs. These arcs are generated by The extent of the hotspot region is clearly visible in these theparticlesinthebrightspot region. Thelower arcisfrom maps. The simulated trailed spectra for two orbital periods the particles in the gas stream that collide with the parti- Spectroscopic simulations of superhumps 7 Figure4. Dissipationmapsfortheaccretiondiscusedintheanalyticalmodel.Thedissipationisdiplayedusingalogarithmicgreyscale forbinary phases 0.00, 0.36, and 0.64respectively. The solidcurve isthe Roche lobeof the primary.The secondary’s orbital motion is vertiallyupwards.Thearrowdisplayedineachmapindicatesafixeddirectiontoanobserver,withphases definedasinFigure1. clesintheaccretiondisc.Theotherarcistheemissionfrom tion in the line profile which occurs in the SPH simulation. theparticlesintheaccretion discthatareinthebrightspot At first glance the pattern appears to repeat from orbit to region. The disc spiral density compression waves produce orbit, but in fact there is a gradual evolution which can be enhanced emission that can be seen in Figure 7 as features clearlyseeninFigures8and9.Thespectralpatternrepeats in the trailed spectrogram along the velocity lines ∼300 after a complete disc precession period. kms−1and ∼−700 kms−1. Figure 8 shows the trailed spec- trogram for a complete disc precession cycle, equivalent to 3.4 Doppler tomography that of thesimple 3D model in Figure 6. Allofoursimulatedtrailedspectrogramsaredominated Dopplertomographywasdevelopedtointerprettrailedspec- bystrongdouble-peakedemission whichexhibitssignificant trograms of CVs (Horne& Marsh 1988; Marsh 2000). Gen- variability.Emissionfromthebrightspotregionisprominent erally observers use a maximum-entropy inversion process in all spectrograms. InFigures 5to 8thesummation of the to generate the most likely Doppler map from the spectra. emission horizontally and vertically is shown, so the upper Dopplertomography condenses theinformation in atrailed panel is the mean line profile, and the right-hand panel is spectrogramintoavelocitymapoftheemissionregionsina the light curve. The light curve in Figure 6 clearly shows system. The observed spectral line profiles as a function of the modulation in the superhump profile which was fully orbital phaseareinverted torevealthebest-fittingDoppler explored in RHP2001 and Rolfe (2001). map of the line flux as a function of velocity in the binary Figure 9is amontage of SPH trailed spectra for differ- frame. It implicitly assumes: ent orbits sampling a complete disc precession. Each spec- • the emission regions are fixed within the co-rotating trogram uses the same linear greyscale and the number in frame the upper left-hand corner of each image indicates the disc • theemissiondoesnotvarywithtimeoverthecourseof precession phase. This figure clearly demonstrates the vari- a binary orbit ability of accretion disc spectra with respect todisc preces- • all points are visible at all times sionphase.WehaveplacedontheInternetaMicrosoft AVI • all motion is in theorbital plane movieshowing thespectrogram’s evolution overa complete • the width of theprofiles from any point is small disc precession cycle (http://physics.open.ac.uk/FHMR/). Inordertocompareourhighresolutiontrailedspectrawith Oursimulations, thosedescribedin (Murray 2000)and currentobservations,wereducedthespectralandtimereso- eclipsemappingofV348Pup(Rolfe et al. 2000)revealdiscs lutionofeachspectrogram ofFigure9byafactorof25and which violate the first two of these assumptions. Hence, theresult is displayed in Figure 10. for precessing non-axisymmetric discs, it is safer to work In all of our trailed spectrograms we see the effects of with trailed spectrograms directly rather than construct- theprecessingdisc.InFigure6weseealongperiodS-wave ingDopplertomograms. Nevertheless,manyof therelevant intheenvelopeofthediscemission,causedbythedisc’spre- observational publications in the last decade have included cession on P . The orbital modulation in the brightspot tomography, so we have constructed tomograms from our prec emission causes the shorter timescale S-wave which is most simulated line profiles. clearly apparent where it crosses each of the two peaks in Doppler maps represent the binary configuration in theline profile. two-dimensional velocity-space. Figure 11 shows an exam- InFigure7weshowthemorecomplex orbitalmodula- ple of an instantaneous Doppler map generated from the 8 S. B. Foulkes et al. Figure 5. Trailed spectrogram using a linear greyscale for the analytical accretion disc model. The spectrogram is for two orbital periods. The curves on the top and right of the spectrogram are the average line flux vertically and horizontally respectively, i.e. they constitute ameanlineprofileandalightcurverespectively. SPHdissipationdata.Theprimaryandsecondaryvelocities cellmultipliedbytheparticle dissipation,isplottedusinga are overlaid on the Doppler map y-axis with the secondary linear grey intensity scale. having a positive y-velocity. The lower of the two overplot- Figure 1(d, e, f) demonstrates that the emission from ted arcs indicates the ballistic trajectory of thegas stream; the disc changes substantially in the course of a binary or- the upper arc indicates the velocity of a circular Keplerian bit, and is far from fixed in the co-rotating frame. Hence disc at positions along the ballistic stream trajectory: the theinstantaneous Dopplermap varies significantly with or- stream’s ‘Kepler shadow’. The velocity of each particle is bital phase. Figure 1(g, h, i) shows detail of instantaneous used to transformed each particle into velocity space. The Dopplermapscorrespondingtothedissipationsnapshotsin velocity density at each velocity cell in the Doppler map Figure 1(d, e, f). is the sum of the dissipation from each particle mapped to Theprominentdissipationneartheletter“G”onFigure thatcell,i.e.thenumberofparticlesoccupyingeachvelocity 1(i) arises from thegas stream as it impacts theouteredge of the accretion disc. This gas is accelerated until it has Spectroscopic simulations of superhumps 9 Figure 6. Trailedspectrogram usingalineargreyscaleforthesimple3Daccretiondiscmodel.Thespectrogramisforacompletedisc precession.Thecurvesonthetopandrightofthespectrogram aretheaveragefluxverticallyandhorizontallyrespectively. the same speed as the gas in the disc, producing the small InFigure1(g)therearethreearcsofenhancedemission. very bright arc close to the impact point. The letter “H” Two of them are produced by the stream particles and the onthisDopplermapindicatesthearcofemission produced disc particles at the impact; the third, uppermost, velocity by the disc particles in the impact region. As the two flows arc is produced by emission from the centrifugally expelled converge (progressing to velocity locations which move in material, which can be seen extending outside the primary ananti-clockwisesense)therelativevelocityofthetwoarcs Roche lobe in Figure 1(d). The material outside the lobe is decreasestozero.Emission from thegasdownstream ofthe moving up and slowly to the right in Figure 1(d), and so it impactpointcanbeseentrailingaway(anti-clockwise)from appears in the upper right quadrant of the velocity map in thebrightspotimpactpoint.Theouterdiscinthesimulation Figure1(g),whereitextendsacrossthedonorvelocitylobe. does not consist of particles executing circular Keplerian The disc in our simulation is eccentric, precessing orbits, so the disc velocity at the impact point does not rapidlywithrespecttotheco-rotatingframe(sinceitispre- generally lie on the Kepler shadow trajectory overplotted cessing only slowly in the inertial frame), and continuously on themap. flexing and relaxing on the superhump period. This causes 10 S. B. Foulkes et al. Figure 7. TrailedspectrogramusingalineargreyscalefortheSPHaccretiondiscmodel.SeeFigure5caption. theappreciable changes intheinstantaneousDopplermaps particles.Thesetwobrightarcsshowthechangeofvelocity produced at different orbital phases, Figure 1(g,h and i). (and position) of the impact region as a function of orbital phase.Thereisalsoanotherlowerintensityarcinthelower Since observers cannot produce an instantaneous left-handquadrant.Thesourceofthisarcisoneofthespiral Dopplermap,weshowinFigure12,middleimage,aDoppler density waves seen in Figure 1(d). dissipationmapforacompleteorbitalperiod.Thismapwas TheSPH trailed spectrogram, shown in Figure12,was generatedbyaveragingtheinstantaneousDopplermapsfor eachmodeltimestep,0.01Ω−1,i.e.629(≈100×2π)instan- inverted using doppler, a maximum entropy Doppler to- orb mographypackagedevelopedbyTomMarsh4,andisshown taneous maps corresponding to a complete orbit were aver- in thebottom image ofFigure 12.Aconstant Dopplermap aged.Twolargeemissionarcsdominatethismap.Thelower wasgeneratedusingmakimgandthenscaledusingoptscl. arcissituated justabovethegasstream ballistictrajectory The Doppler map had a starting chi-squared value of 26.6 and was generated by the gas stream impacting the outer edge of the accretion disc. The merging streams generate a large dissipation. The upper arc is generated by the parti- cles in the accretion disc that interact with the gas stream 4 http://www.astro.soton.ac.uk/∼trm/tar/doppler.tar.gz