双星观测量计算方法的改进及应用∗

陈 健夏 芳 傅燕宁

(1中国科学院紫金山天文台南京210034)

(2中国科学院大学北京100049)

1 引言

近几十年来,通过不断的观测和拟合研究,发现并确认的恒星双星系统越来越多[1].目前认为有超过一半的恒星是双星[2],且其中大部分系统都有亮子星,但由于观测条件的限制,目前已知轨迹的系统只有约7000个[3],在已经发现的双星中占的比例不到5%.随着观测技术的发展,高精度的恒星运动学观测资料大量积累.华盛顿双星星表(The Washington Double Star Catalogue,WDS)中已经收录了142645个双星系统的1619023个相对位置测量数据1http://ad.usno.navy.mil/wds/wds.html,这使得双星轨道拟合工作发展迅速,其给出的大批量的系统亮子星的位置,丰富了亮星星表参考架的内容[1].同时,双星轨道拟合给出的恒星动力学质量参数是恒星系统动力学研究和恒星物理演化研究必不可少的参量,也为恒星经验质光关系的研究提供了更多精确的参考样本[4−6].

双星系统的基本动力学模型是二体问题,其对应的运动方程是可积的,积分常数一般取为质心的位置速度和二体相对运动的轨道参数[1],在双星各自的质量已知的情况下,就可以预测双星在任意一个时刻的位置和速度.反之,在已知多个时刻的观测位置和速度之后,可以通过拟合得到积分常数和质量.针对双星的观测资料主要有相对位置资料、视向速度资料以及依巴谷、Gaia等天体测量数据.尽管二体运动学模型具有足够的完备性,但此前在针对这些观测资料的拟合工作中,计算观测量理论值时唯一考虑的改正是由岁差引起的切平面正北方向的改正[7].而目前有研究者们注意到现有观测量计算模型所依赖的切平面局部参考系本身就存在一些近似问题[3,8].首先,双星轨道的空间指向由切平面给出[9],切平面在精度较低时可以粗略地定义为垂直于某颗子恒星的视线方向,但我们知道,这个视线方向由于受到恒星自行的影响不断发生变化,因此,在高精度恒星定位研究中需要给切平面一个精确的定义或者采用其他更精确的方法描述双星轨道空间指向;其次,由于空间透视效应的存在[8],轨道运动在切平面上的平行投影并不是我们真正看到的图像,而是假想的从无穷远处观测双星所看到的图像.Kaplan[3]2015年在IAU大会指出这些问题导致的相对位置计算偏差有可能达到40 mas.早期因为相对位置测量误差普遍大于0.1′′,所以近似的问题一直没有引起重视.而目前地面的干涉观测仪器的精度已经达到或好于1 mas[1],为避免近似带来的误差,采用严格的立体几何方法计算双星的观测量是必须的.

本文第2节介绍观测量计算方法改进的具体内容;第3节介绍改进后的方法在双星观测量计算与轨道拟合上的应用;第4节讨论观测量计算方法改进的意义及未来继续提高模型精度需要考虑的其他问题.

2 改进观测量计算方法

我们以牛顿二体模型为基础,通过严格的立体几何推导建立从模型参数到观测量的计算方法.

2.1 模型参数的选取

α、δ、ϖ分别为双星质心在t0时刻的赤经、赤纬、视差.µα∗、µδ分别为双星质心在t0时刻自行的赤经和赤纬分量.vr0为双星质心在t0时刻的视向速度.这些量都是描述双星质心运动的参数.

a0是双星轨道的半长径,以au为单位.e是轨道偏心率.τ是轨道上天体过近心点的时刻.P是轨道周期.这4个参数用于描述二体在轨道平面上的相对运动.q是双星的主星质量mA与伴星质量mB之比,用于将相对运动归算到质心.以上5个参数都不会因采用坐标系的不同而变化.需要注意的是,参数a0在旧模型中通常采用的定义是a=a0ϖ,为了保持一致,新模型的结果最后会由a0换算到a.

最后介绍用于确定双星轨道指向的3个参数:

(1)轨道倾角i0:轨道法向与ICRS北极方向所成的角.

(2)升交点经度Ω0:将ICRS原点平移至主星,Ω0就是该参考系中二体相对运动轨道关于赤道面升交点的赤经.

(3)近心点角距ω0:升交点沿二体轨道顺行到近心点所走过的角度.

在近似的观测量计算方法中上述3个指向角(分别记为i、Ω、ω)是关于“切平面”给出的.其中忽略了切平面的时变性.Kaplan讨论了这种近似对子星运动状态计算的影响[3],本文则以子星相对位置和视向速度这两种常见的观测量为例,进一步讨论其对轨道拟合的影响.

2.2 严格计算观测量

下面给出从模型参数到任一时刻t的子星相对位置和视向速度的计算过程.

首先由α、δ、ϖ、µα∗、µδ、vr0计算出t0时刻双星质心相对太阳系质心的位置速度,再根据质心运动定理求出t时刻双星质心相对太阳系质心的位置rC(t)和速度(t).

接着计算轨道位置,由模型参数a0、e、τ、P求出t时刻双星在轨道坐标系中的位置及速度,再由i0、Ω0、ω03个角度旋转至ICRS的坐标方向下的相对位置速度.然后根据质量比q归算到相对质心的位置速度.

综合上面的结果,求出两颗星在t时刻分别相对于太阳系质心的位置rA(t),rB(t)和速度(t),(t).至此,子星的运动状态就完全确定了,由rA(t),rB(t),(t),(t)即可得到任意一种运动学观测量.

相对位置(ρ(t),θ(t))是双星观测中常用的观测量,其中ρ就是rA(t),rB(t)两个位矢的夹角:

之所以(1)式中通过求rA(t),rB(t)的外积而非内积来解ρ,是因为ρ的值非常小(1′′量级),这种情况下sin ρ比cos ρ对ρ的变化更敏感.

θ是天球上AB的方向关于AP的方向的方位角,如图1所示,O为坐标原点,θ在3维空间中实际上是平面ABO与平面APO所成的角.为了得到θ与rA(t),rB(t)的关系,我们在图1中引入了几个辅助量,其中为ICRS的z轴方向为平面ABO的法方向,为平面APO的法方向.易得∠GOF= θ,又由和可得θ满足的关系式:

其中ez={0,0,1}是北极方向的单位矢量.可以从图1中看出θ就是大圆GFED上F关于G的方位角,而因为在ICRS赤道面上,则通过在z轴上分量的正负可以判断出θ小于或大于180◦.

除了相对位置(ρ(t),θ(t)),观测量还有两子星各自的视向速度vrA(t)和vrB(t),它们就是两子星速度在其各自的位置rA(t),rB(t)方向上的分量:

图1 θ在天球上的示意图.其中O为坐标原点,即太阳系质心,A、B分别为主星、伴星在天球上的位置,P为北天极,大圆EF GD以A为极Fig.1 The diagram of θ on the celestial sphere.O is the origin.A and B are the primary and secondary components on the celestial sphere,respectively.P is the north pole.A is the pole of great circle EF GD

最后根据以上过程求观测量对每一个模型参数的一阶偏导数,再用BVLS(Bounded Variable Least Squares algorithm)[10]方法拟合.需要注意的是,相对位置和视向速度对质心参数α、δ、µα∗、µδ的偏导数都很小,即相关性很弱,所以这4个参数不参与拟合,它们应参考目前精度较高的观测资料给定.这种严格的观测量计算方法涉及到全局量与局部量的加减,因此存在有效数字位数损失的问题,我们在程序中对所有浮点数采用4倍精度,以避免精度损失对结果的影响.

2.3 严格的观测量计算方法的优势

严格的观测量计算方法的优势在于从全局参考系出发描述双星的运动,解决了近似的观测量计算方法中回避了的透视问题及坐标架指向不固定问题.图2显示出切平面局部参考系未考虑透视所引起的偏差,其中B′为伴星在切平面上的近似投影位置,B′′为伴星在切平面上的实际投影位置.切平面指向不固定主要由双星系统的自行导致,引起的观测量偏差通常会随着时间的推移而逐渐增大,必须通过质心运动才能改正.

图2 透视效应对双星投影位置的影响.其中A、B分别为主星和伴星,B′B平行于A的视线方向且交切平面于点B′,B′′为B的视线方向与切平面的交点Fig.2 The prospective effect on the projection of binary.In this graph,A and B are the primary and secondary components,respectively.B′B is parallel to the line of sight to A,and intersects the plane of sky in B′,while the line of sight of B intersects the plane of sky in B′′

从近似的观测量计算方法产生偏差的机制来看,双星离我们观察者越近,并且双星间的距离越大,产生的偏差也就越大,因此观测量计算方法的改进对于太阳系邻近区域中的远距双星更有帮助.

3 严格的观测量计算方法的应用与结果分析

上一节我们介绍了严格的观测量计算方法,本节主要讨论其在实际系统中的应用.

首先,我们在第6目视双星星表(Sixth Catalog of Orbits of Visual Binary Stars)中找出所有含亮子星的双星系统,并统计严格的与近似的观测量计算方法计算出的相对位置的差异,发现有不少系统的差异已经超过了当前的观测误差,其中差异较大的系统大多具有较大的轨道半长径和视差.

从差异较大的系统中我们选取了观测资料相对丰富的半人马座α双星作为拟合的试验对象.该系统是离太阳系最近的亮双星,自行与半长径均较大,因此采用近似观测量计算方法造成的偏差也较大.针对它的轨道拟合工作非常多[11–13],最新的轨道参数是Kervella等人给出的,所用观测数据为1752–2016年的相对位置数据和1900年以来的视向速度数据[13],我们用了和他们同样的观测数据,拟合中的权重也依据Kervella等人文中给出的误差来分配.拟合中需要固定不动的α、δ、µα∗、µδ4个质心参数中,µα∗、µδ也参考了Kervella等人的工作[13],α、δ则是由依巴谷星表中子星(HIP71683和HIP71681)各自的位置计算得来.接下来通过BVLS方法分别用严格的与近似的观测量计算方法进行拟合.得到的χ2非常接近,但是通过两种方法解出的参数不一样(见表1),这说明两种观测量计算方法虽然与观测的符合程度相似,但调整参数的方向却不同.近似的观测量计算方法可能使模型参数向着错误的方向调整了.于是我们进一步对参数的解进行分析.

参考实际观测数据的精度和轨道覆盖率,我们产生了100组仿真数据,并对这些数据进行再次拟合,最后在表1中给出参数误差,并列出了误差范围内出现真值的可能性,作为参数准确性的判断依据.从表1中我们可以看出,依靠严格的观测量计算方法解出的参数的准确性均高于近似方法.尤其对vr0和q两个参数,采用近似的观测量计算方法拟合,能恢复到真值的机率只有不到1%,这是因为这两个参数与观测量vrA(t)、vrB(t)的相关性很强,依靠严格的观测量计算方法可以把vr0和q定得很精确,而近似的观测量计算方法在描述双星运动方面是有偏差的,参数解与真值实际上差别很大.

我们在图3中显示了用两种观测量计算方法拟合仿真数据所得参数结果的相对偏差量统计.横坐标上每一项代表一个参数的拟合值与真值的相对偏差量,对应两组统计图,其中空心的是严格的观测量计算方法拟合值的偏差量统计,实心的代表近似方法.每一组统计图中,箱内的一根横线表示这组相对偏差的中位数,小方块表示平均数,箱代表这组相对偏差样本按大小排序后位于中间的50%所在的范围,误差棒则代表位于中间的98%的样本所在的范围,可以看出,中位数与平均数基本都在箱图的中间,它们与零点的差距越大说明对应的拟合结果越不合理.从图3中可以看出,对每一个参数,严格的观测量计算方法拟合结果总是比近似的更靠近零点,即更合理.

表1 半人马座α双星的参数解及其准确性Table 1 The fitted parameters of α Centauri and their confidence levels

图3 各参数拟合结果偏差图.空心框图与实心框图分别对应严格的与近似的观测量计算方法的仿真拟合统计结果Fig.3 The relative deviations of parameters.The open boxes and the filled boxes,respectively,show the statistical results of fitting simulated data based on the exact and approximate algorithms for computing the observables

图3中看出近似的观测量计算方法拟合参数解的相对偏差最大不超过1%,但这并不代表参数解的偏差可以被忽略.我们用表1中两种方法拟合的参数分别计算4种观测量并求差值,将这个差值的变化曲线画在图4中,可以看出这个偏差将会导致预报的双星相对位置偏差达到15 mas,视向速度偏差达到12 m·s−1,这都远超过了当前这两种观测数据的观测误差.同时,随着时间的推移这些偏差都有逐步增大的趋势.

图4 严格的与近似的观测量计算方法对半人马座α观测数据的拟合结果在计算观测量上的偏差.其中ρ、θ为相对位置的坐标分量,vrA、vrB为主星和伴星的视向速度,下标为“new”的是严格的观测量计算方法拟合出的参数计算出的观测量,下标为“old”的是近似方法拟合出的参数计算出的观测量Fig.4 The difference between observed quantities computed from parameters fitted by the exact and approximate algorithms to the observational data of α Centauri. ρ and θ are the components of the relative positions.vrA and vrB are the radial velocities of the major component star and the minor one,respectively.The subscript“new”means that the observables are computed from the parameters fitted by the exact algorithm,and the subscript“old” means that the observables are computed from the parameters fitted by the approximate algorithm

4 结论和展望

随着观测数据精度的不断提高,采用严格的观测量计算方法是必须的.本文对此前观测量计算方法中忽略的两个近似因素进行了分析,并给出了改进后的观测量计算方法.实际系统的应用结果表明,尽管严格与近似的观测量计算方法与观测数据的符合程度相当,但近似方法得到的拟合参数偏离真值的程度更大.尤其是视向速度和质量比这两个在恒星动力学研究和物理研究方面都很重要的参数,近似方法给出正确值的可能性很低.另外,这种偏离在恒星观测量计算中体现的偏差也不可忽视.因此,对双星系统,尤其是太阳系附近的远距双星系统,需要采用严格的观测量计算方法来拟合.

当然,关于观测量计算模型的改进,今后还有一些其他因素需要考虑.首先是光行差引起的双星相对位置的偏差,其大小约为两子星角距离的10−4[3],即角距离10′′的双星偏差为1 mas.其次是子星离观测者距离不同导致光线传播的时间不同,观测到的双星位置对应不到同一时刻,该因素对轨道倾角大、空间距离大的双星影响较大[8].综上所述,在提高双星运动模型精度方面还有大量的工作需要完成.

模型精度的提高需要大量高精度观测数据的支持,观测技术的不断发展为之提供了可能.目前最具代表性的就是欧洲空间局(ESA)于2013年发射的Gaia卫星,它对亮于15 mag的恒星能提供精度达到7–20µas的天体测量数据[14].除了这些高精度的新观测数据,长期的历史观测资料在双星轨道拟合方面的作用也不容忽视[15].年代久远的天文观测底片也可以通过扫描技术重新得到整理,使位置测量精度提高到0.1′′左右[16].有了高精度的模型和数据,我们可以给出大批量双星的新参数.在此基础上,可以进一步开展的工作有:双星系统轨道性质的统计研究;编制和维持更高精度、更高网格密度的星表;利用各种轨道拟合过程得到恒星质量等物理参数;研究多星系统和星团的稳定性或者动力学演化过程;检验恒星的演化模型等[1].

致谢 感谢Pierre Kervella与Dimitri Pourbaix为我们提供半人马座α的观测资料.