التطور طويل الأمد للأنظمة الكوكبية ذات كوكب أرضي وكوكب عملاق

Nikolaos Georgakarakos,1 Ian Dobbs-Dixon,1 Michael J. Way2
1New York University Abu Dhabi, Saadiyat Island, P.O. Box 129188, Abu Dhabi, UAE
2 NASA, Goddard Institute for Space Studies, New York, USA
E-mail: ng53@nyu.edu, georgakarakos@hotmail.com
(قُبل XXX. استُلم YYY؛ بصيغته الأصلية ZZZ)
الملخص

ندرس التطور المداري طويل الأمد لكوكب أرضي تحت تأثير الاضطرابات الثقالية الناجمة عن كوكب عملاق. وعلى وجه الخصوص، نهتم بالحالات التي يكون فيها الكوكبان في المستوى نفسه ومتقاربين نسبيا. ونفحص التهيئتين الممكنتين: أن يكون مدار الكوكب العملاق إما خارج مدار الكوكب الأصغر أو داخله. ويُوسَّع كمون الاضطراب إلى رتب عالية، ثم يُشتق حل تحليلي لمدار الكوكب الأرضي. بعد ذلك تُقارَن التقديرات التحليلية بنتائج التكامل العددي لمعادلات الحركة الكاملة، ونجد أن الحل التحليلي يعمل على نحو معقول. ومن النتائج اللافتة أن التقديرات التحليلية الجديدة تحسن كثيرا تنبؤات المقاييس الزمنية للتطور المداري للكوكب الأرضي مقارنة بتوسع من رتبة الثماني الأقطاب. وأخيرا، نناقش بإيجاز التطبيقات المحتملة للتقديرات التحليلية في مسائل الفيزياء الفلكية.

keywords:
الميكانيكا السماوية - الكواكب والأقمار: التطور الديناميكي والاستقرار – الكواكب والأقمار: الكواكب الأرضية
pubyear: 2015pagerange: التطور طويل الأمد للأنظمة الكوكبية ذات كوكب أرضي وكوكب عملاق7

1 المقدمة

قبل ثلاثة عقود، كانت الكواكب الوحيدة التي نعرفها هي تلك الموجودة في نظامنا الشمسي. أما اليوم، فقد اكتُشف أكثر من 3400 كوكب في 2551 نظاما كوكبيا، مع وجود 581 من تلك الأنظمة تحتوي على كوكبين على الأقل (www.exoplanet.eu). وعلى خلاف نظامنا الشمسي، حيث تمتلك الكواكب، باستثناء عطارد، مدارات شبه دائرية ومشتركة المستوى، فإن الأنظمة الكوكبية الخارجية التي اكتُشفت حتى الآن تمتلك اختلافات مركزية كوكبية تمتد من الصفر إلى ما يقارب الواحد، كما في حالة HD20782b (O’Toole et al., 2009). وقد تكون لدى أنظمة أخرى مدارات مائلة بعضها بالنسبة إلى بعض بوصفها تهيئة ممكنة، كما يبين اكتشاف υ Andromedae حيث يمتلك الكوكبان العملاقان ميلا متبادلا قدره 30 (McArthur et al., 2010). وللاطلاع على مراجعة حديثة جيدة عن بنية الأنظمة الكوكبية الخارجية، يمكن الرجوع إلى Winn & Fabrycky (2015).

عادة ما يكون تقريب الثلاثي الهرمي نهجا جيدا لنمذجة السلوك طويل الأمد لنظام ثلاثي الأجسام، إذ يمكن النظر إلى النظام على أنه ثنائيان متفاعلان؛ فالجسمان المتقاربان يشكلان الثنائي الداخلي، بينما يكون الجسم الثالث على مدار أوسع بالنسبة إلى مركز كتلة الثنائي الداخلي. في مثل هذه الحالة، يمكن التعبير عن كمون الاضطراب على هيئة متسلسلة قوى في نسبة نصفي المحورين الرئيسيين، ويمكن قطعها عند أي رتبة تُعد ضرورية للمسألة المدروسة. وفي الأنظمة النجمية أو الكواكب الموجودة في ثنائيات، من المعتاد تضمين الحدود الناشئة من الحدين P2 وP3 عند توسيع الهاملتوني الاضطرابي بدلالة كثيرات حدود لوجندر، وهي ما تسمى بحدي رباعي القطب وثماني القطب؛ انظر مثلا Ford et al. (2000)، Lee & Peale (2003)، Georgakarakos (2004)، Naoz et al. (2013). ويكون هذا النهج عادة كافيا عندما يوجد نجمان على الأقل في نظام الأجسام الثلاثة، أو كوكبان في مدارات عالية الميل. ويعود ذلك إلى أن الفواصل يجب أن تكون كبيرة بما يكفي ليبقى النظام مستقرا (انظر مثلا Georgakarakos, 2013). غير أنه في نظام شبه مشترك المستوى يحوي نجما واحدا فقط، يمكن للكواكب أن تكون أكثر تقاربا بكثير وأن تبقى مع ذلك مستقرة.

في سلسلة من الأوراق (Georgakarakos, 2002, 2003, 2009) حصلنا على تعبيرات تحليلية للاختلاف المركزي وخط طول الحضيض في الأنظمة الثلاثية الهرمية غير الرنينية والمشتركة المستوى. وقد أُجريت تلك الحسابات على مقاييس زمنية قصيرة (من رتبة الفترات المدارية) وعلى مقاييس زمنية دهرية. وعلى الرغم من عدم وجود أي قيد صريح على نسب الكتل في الأنظمة، كان ذلك العمل موجها إلى التركيز على الأنظمة النجمية. وفي Georgakarakos (2006) اختبرنا صلاحية نتائج Georgakarakos (2002, 2003) للأنظمة الكوكبية، ولكن من حيث نسب الكتل فقط، وصولا إلى أجسام بحجم المشتري فحسب. في هذا العمل نحاول توسيع عملنا السابق باشتقاق تعبيرات تحليلية لحركة نظام مؤلف من كوكبين تقع جميع أجسامه في المستوى نفسه. ونولي اهتماما خاصا للحالة التي يكون فيها أحد الكوكبين أصغر كثيرا من الآخر، لأن هذه التهيئة يمكن أن تكون لها تطبيقات كثيرة في دراسات قابلية السكن. ومن ثم فإن اهتمامنا الأساسي في العمل الحالي ينصب على أنظمة تحوي كوكبا شبيها بالأرض أو أرضا فائقة وكوكبا عملاقا. وقد تكون النتائج قابلة للتطبيق أيضا على حركة أنظمة أخرى، مثل نظام يضم نجما وكوكبا عملاقا وقمره الخارجي. وأخيرا، نود الإشارة إلى أننا في العمل الحالي سنركز فقط على الحركة طويلة الأمد للنظام.

تأتي بنية الورقة على النحو الآتي: في القسم 2، نعرض صياغة المسألة والحل التحليلي لها. وفي القسم 3، نختبر النتائج الجديدة باستخدام تكاملات عددية لمعادلات الحركة الكاملة، وفي القسم 4 نلخص نتائجنا ونناقش الآفاق المستقبلية.

2 النظرية

كما ذكرنا بإيجاز في المقدمة، يهدف هذا البحث إلى اشتقاق تقديرات تحليلية لحركة كوكب بكتلة الأرض أو الأرض الفائقة تحت التأثير الثقالي لكوكب عملاق، مع وقوع جميع الأجسام في المستوى نفسه. وسنركز على التطور المداري طويل الأمد، لأن مساهمة الحدود قصيرة الفترة في التطور المداري لنوع الأنظمة التي ندرسها صغيرة جدا (باستثناء مدارات الكوكب العملاق شبه الدائرية). وبما أن أنصاف المحاور الرئيسية في المسألة الدهرية لا تتغير (انظر مثلا Harrington, 1968)، فإننا لا نحتاج إلا إلى تحديد الاختلافات المركزية وخطوط طول الحضيض لكواكبنا كي نحصل على وصف ديناميكي كامل لنظامنا. وبافتراض أن مدار الكوكب الأكبر لا يتأثر بوجود الكوكب الأصغر، يكفي أن نعالج التطور المداري للأخير.

نبحث إعدادين مختلفين. في الأول، يتكون الثنائي الداخلي في ثلاثينا الهرمي من النجم والكوكب الأصغر، بينما يكون الجسم الأكبر في مدار خارجي بالنسبة إلى ذلك. وسنسمّي تلك الحالة فيما يأتي الحالة الداخلية. وفي الثانية المسماة الحالة الخارجية، ندرس الوضع الذي يكون فيه مدار الكوكب العملاق داخل مدار الكوكب الأصغر. وهذه أيضا تهيئة ممكنة كما تبينها دراسات حديثة (Raymond et al., 2006; Ogihara et al., 2014) واكتشاف أنظمة مثل Kepler-87 (e.g. see Ofir et al., 2014).

فيما يلي، سنستخدم الصياغة الهاملتونية لدراسة الحركة الدهرية للأنظمة قيد البحث. أولا، سنشتق الهاملتوني الدهري ثم معادلات الحركة التي تُحل بعد ذلك تحليليا، مانحة إيانا حلا لمسألتنا. ويُفترض أن مدار الكوكب الأصغر دائري في البداية وأن اختلافه المركزي الأعظمي سيبقى منخفضا إلى حد ما. إضافة إلى ذلك، في الحل المداري للحالة الخارجية، نحتفظ بحدود من الرتبة الثامنة في اختلاف مركزية الكوكب العملاق، لأن معظم المشتريات الحارة أو الدافئة لها اختلافات مركزية صغيرة إلى متوسطة. وتُعامل جميع الأجسام ككتل نقطية، وفي هذه المرحلة لا ننظر إلا في التفاعلات الثقالية على المستوى النيوتني. وبالنسبة إلى جميع العناصر المدارية، يشير الدليلان 1 و2 إلى الكوكب الداخلي والخارجي على الترتيب.

2.1 توسيط الهاملتوني

في عملنا السابق (Georgakarakos, 2002, 2003)، أُنجز اشتقاق الهاملتوني الدهري بطريقة فون زيبل، التي تستخدم تحويلات قانونية لإزالة المتغيرات السريعة. وإلى جانب إعطاء هاملتوني كامل من حيث الاختلافات المركزية، تنتج طريقة فون زيبل أيضا حدودا من الرتبة الثانية في الكتل، وهو أمر لا يمكن بلوغه بطرائق أخرى كثيرة (مثل طريقة «المقصات»). غير أن هذه السمة الأخيرة لطريقة فون زيبل مهمة في الأنظمة ذات الكتلة الخارجية الكبيرة، وهو ليس الحال هنا. ومن ثم، ولغرض هذه الورقة، نعتمد طريقة مختلفة للحصول على الهاملتوني الدهري.

اتباعا لـ Lee & Peale (2003)، نستخدم المتغيرات الآتية، وهي مجموعة من المتغيرات القانونية (Brouwer & Clemence, 1961; Murray & Dermott, 1999)، وتُعرّف لنظام ثلاثي الأجسام مشترك المستوى كما يأتي:

L1=mn1a12,l1,L2=n2a22,l2,G1=L11e12,ϖ1,G2=L21e22,ϖ2, (1)

حيث إن ai وei وϖi وli هي أنصاف المحاور الرئيسية، والاختلافات المركزية، وخطوط طول الحضيض، والشذوذات المتوسطة للمدار الداخلي (i=1) والمدار الخارجي (i=2) على الترتيب. علاوة على ذلك،

m=m0m1m0+m1and=m2(m0+m1)M,

هي ما يسمى بالكتل المختزلة، حيث إن m0,m1 وm2 هما كتلتا النجم والكوكبين على الترتيب. والكتلة الكلية للنظام هي M=i=02mi. وباستخدام هذه المتغيرات، يُكتب الهاملتوني لنظام ثلاثي هرمي على الصورة

H=H0+H1+Hp, (2)

حيث

H0=𝒢2m03m132(m0+m1)L12 (3)

هي الطاقة الكبلرية للمدار الداخلي،

H1=𝒢2(m0+m1)3m232ML22 (4)

هي الطاقة الكبلرية للمدار الخارجي، و

Hp=𝒢m2(m0+m1Rm0r02m1r12), (5)

هو الهاملتوني الاضطرابي، حيث إن r02 وr12 هما المسافتان بين m0 وm2 وبين m1 وm2 على الترتيب. و𝒢 هو ثابت الجاذبية، وR هو المتجه الذي يصل مركز كتلة الثنائي الداخلي بالجسم الثالث، بينما سيشير r إلى متجه الموضع النسبي للمدار الداخلي (انظر الشكل 1).

Refer to caption
Figure 1: تحليل جاكوبي لمسألة الأجسام الثلاثة. لاحظ أن الحجم النسبي لكتل الكواكب يدل على الحالة الخارجية، لكن التهيئة النسبية صالحة لكلتا الحالتين.

يمكن كتابة المعادلة (5) على الصورة

Hp=𝒢m2(m0+m1Rm0|R+μ1r|m1|Rμ0r|), (6)

مع

μi=mim0+m1,i=0,1.

وبما أننا نتعامل مع نظام هرمي وr/R<<1، يمكن تمثيل مقلوبات المسافات بدلالة كثيرات حدود لوجندر Pn كما يأتي (انظر مثلا Murray & Dermott, 1999):

m0|R+μ1r|=m0Rn=0(μ1rR)nPn(cosθ) (7)

و

m1|Rμ0r|=m1Rn=0(μ0rR)nPn(cosθ), (8)

حيث إن θ هي الزاوية بين المتجهين r وR. بعد ذلك، وبعد بعض المعالجات الجبرية الأساسية، تصبح المعادلة (6):

Hp=𝒢m0m1m2Rj=2Mj(rR)j𝒫j(cosθ), (9)

حيث

Mj=m0j1(m1)j1(m0+m1)j.

لوصف الحركة طويلة الأمد للنظام، يجب توسيط الهاملتوني أعلاه على الزوايا السريعة، بعد إدراج التعبيرات المناسبة لمتجهي الموضع النسبيين والزاوية θ. وتُعطى كثيرات حدود لوجندر بـ

Pn(cosθ)=2nk=0ncoskθ(nk)(n+k12n), (10)

أما بقية الكميات فنستخدم لها التعبيرات الآتية:

rR = a1(1e1cosE1)1+e2cosf2a2(1e22) (11)
cosθ = cos(f1+ϖ1f2ϖ2) (12)
cosf1 = cosE1e11e1cosE1 (13)
sinf1 = 1e12sinE11e1cosE1, (14)

حيث إن E1 هو الشذوذ اللامركزي للمدار الداخلي، وfi يدل على الشذوذ الحقيقي (i=1,2). ومن أجل إيجاد الهاملتوني الدهري، نحتاج إلى ترشيح التأثيرات ذات المقاييس الزمنية من رتبة الفترتين الكوكبيتين. وتمثل هذه التأثيرات في هاملتونيّنا بالشذوذ اللامركزي الداخلي E1 والشذوذ الحقيقي الخارجي f2؛ ومن ثم يمكن الحصول على الهاملتوني الدهري بحساب التكامل المزدوج الآتي:

Hl=H0l+H1l+14π202π02π(1e1cosE1)(1e22)3/2(1+e2cosf2)2HpdE1df2, (15)

حيث يدل الدليل l على التطور طويل الأمد. وعلى الرغم من أننا، لأسباب تتعلق بالتبسيط، لن نستخدم الدليل l فيما يأتي، فإن الهاملتوني والعناصر المدارية سيشيران إلى الكميات طويلة الأمد.

ضمن هذا الإطار، يمكن اشتقاق الهاملتوني المتوسط حتى أي رتبة نختارها. وقد اخترنا في هذا العمل حساب الهاملتوني حتى الرتبة الحادية عشرة. هدفنا الرئيس هو وصف حركة كوكب أرضي تحت التأثير الثقالي لكوكب عملاق عندما تكون نسبة نصفي المحورين الرئيسيين صغيرة نسبيا (a2/a1 قد تبلغ نحو  1.6، 1.7). وقد أظهرت الاختبارات العددية، حتى الرتبة الخامسة عشرة، فروقا مهملة مقارنة بتوسع الرتبة الحادية عشرة لهذا النوع من نسب أنصاف المحاور الرئيسية، مما يشير إلى أن الأخير يقدم تقديرا تحليليا جيدا بما يكفي لنسب أنصاف المحاور الرئيسية والكتل المدروسة في هذا العمل. ومع ذلك، تُعرض التقديرات التحليلية بطريقة تتيح قطعها عند أي رتبة دون الحادية عشرة. ويُعطى مزيد من التبرير لهذا الاختيار في القسم 3.

بما أن مدار الكوكب الأصغر دائري في البداية ويُفترض أنه لن يصبح عظيم الاختلاف المركزي، فإننا، بدلا من اشتقاق الهاملتوني الكامل ثم إهمال حدود O(ei2) في معادلات الحركة الدهرية كما فُعل في Georgakarakos (2003)، نختار إجراء تقريب مكافئ باشتقاق هاملتوني دهري تقريبي تُهمل فيه الحدود من الرتبة O(ei3) وما فوقها. وستقود معادلات الحركة المبنية على الهاملتوني التقريبي في النهاية إلى الحل نفسه الذي نحصل عليه لو أزلنا حدود الرتبة الثانية من معادلات الحركة. غير أننا بهذه الطريقة نتجنب التعبيرات المطولة للهاملتوني ومعادلات الحركة اللاحقة، مع الحفاظ على دقة تقديراتنا التحليلية. ويمكن العثور على الهاملتوني الدهري الاضطرابي الكامل (في الاختلافات المركزية) من الرتبة الحادية عشرة، الذي سيُستخدم لأغراض الاختبار في القسم 3، في الملحق B، وهو متوافق جيدا مع النتائج التي حصل عليها مؤلفون آخرون (مثلا Migaszewski & Goździewski, 2008; Laskar & Boué, 2010).

2.2 الحالة الداخلية

نبدأ باشتقاق الحل التحليلي للحالة الداخلية، حيث يكون مدار الكوكب الأصغر داخل مدار الجسم الكوكبي الأكبر. ومن المعادلة (15)، ومع التعبير بالعناصر المدارية حتى الرتبة الحادية عشرة، يصبح الهاملتوني الاضطرابي طويل الأمد التقريبي:

Hp=Hp2+Hp3+Hp4+Hp5+Hp6+Hp7+Hp8+Hp9+Hp10+Hp11, (16)

مع

Hp2 = 14𝒢m0m1m2(m0+m1)(1e22)3/2a12a23(1+32e12) (17)
Hp3 = 1516𝒢m0m1m2(m0m1)(m0+m1)2(1e22)5/2a13a24e1e2cos(ϖ1ϖ2) (18)
Hp4 = 964𝒢m0m1m2(m03+m13)(m0+m1)4(1e22)7/2a14a25[1+32e22+(5+152e22)e12+354e12e22cos(2ϖ12ϖ2)] (19)
Hp5 = 10564𝒢m0m1m2(m04m14)(m0+m1)5(1e22)9/2a15a26(1+34e22)e1e2cos(ϖ1ϖ2) (20)
Hp6 = 25256𝒢m0m1m2(m05+m15)(m0+m1)6(1e22)11/2a16a27[1+5e22+158e24+(212+1052e22+31516e24)e12+(1894+1898e22)e12e22× (21)
cos(2ϖ12ϖ2)]
Hp7 = 47252048𝒢m0m1m2(m06m16)(m0+m1)7(1e22)13/2a17a28(1+52e22+58e24)e1e2cos(ϖ1ϖ2) (22)
Hp8 = 12258192𝒢m0m1m2(m07+m17)(m0+m1)8(1e22)15/2a18a29[12+214e22+10516e24+3532e26+(9+1892e22+9458e24+31516e26)e12+ (23)
(2974+4954e22+148564e24)e12e22cos(2ϖ12ϖ2)]
Hp9 = 242558192𝒢m0m1m2(m08m18)(m0+m1)9(1e22)17/2a19a210(1+214e22+358e24+3564e26)e1e2cos(ϖ1ϖ2) (24)
Hp10 = 396932768𝒢m0m1m2(m09+m19)(m0+m1)10(1e22)19/2a110a211[12+9e22+1898e24+1058e26+315256e28+(554+4952e22+1039516e24+ (25)
+577516e26+17325512e28)e12+(7154+50058e22+2502564e24+5005128e26)e12e22cos(2ϖ12ϖ2)]
Hp11 = 455262144𝒢m0m1m2(m010m110)(m0+m1)11(1e22)21/2a111a212(2079+18711e22+1309774e24+21829516e26+130977128e28)e1e2× (26)
×cos(ϖ1ϖ2)

بعد إيجاد الهاملتوني، نستطيع الآن اشتقاق معادلات الحركة لنظامنا الثلاثي الأجسام. وباستخدام مجموعتنا من المتغيرات القانونية، تأخذ معادلات هاملتون الشكل الآتي لنظام مشترك المستوى:

dLidt=Hli,dlidt=HLidGidt=Hϖi,dϖidt=HGi,i=1,2. (27)

وبما أن العناصر المدارية المعتادة ليست مجموعة من المتغيرات القانونية، يمكننا استخدام العلاقات أعلاه للحصول على معادلات الحركة بدلالة العناصر المدارية (تفاصيل أكثر في الملحق A). ومع تذكر أننا افترضنا أن مدار الكوكب الأكبر لا يتأثر بوجود الكوكب الأصغر، وأنه لا يوجد تغير دهري في أنصاف المحاور الرئيسية (لأن الهاملتوني مستقل عن li)، فلن نحتاج إلا إلى تعبيرات للاختلاف المركزي وخط طول الحضيض لمدار الكوكب الأصغر. وهذه هي:

de1dt = (m0+m1)1/2(1e12)1/2𝒢1/2m0m1a11/2e1Hϖ1 (28)
dϖ1dt = (m0+m1)1/2(1e12)1/2𝒢1/2m0m1a11/2e1He1. (29)

غير أنه، بدلا من استخدام e1 وϖ1، من الأفضل استخدام المتغيرين x1=e1cosϖ1 وy1=e1sinϖ1، أي مركبتي متجه الاختلاف المركزي الداخلي طويل الأمد. وبهذه الطريقة يمكننا تجنب التفردات في معادلات الحركة عندما يكون الاختلاف المركزي للجسم الأصغر صفرا في البداية. وباستخدام تعبيرات الهاملتوني طويل الأمد التي اشتققناها سابقا، نحصل على نظام المعادلات التفاضلية الآتي:

x˙1 = Asin2ϖ2x1+(B2Acos2ϖ2)y1+Csinϖ2 (30)
y˙1 = (B2Asin2ϖ2)x1Asin2ϖ2y1Ccosϖ2, (31)

حيث إن A,B وC ثوابت وتُعطيان كما يأتي:

A=Ap4+Ap6+Ap8+Ap10 (32)
Ap4 = 315128Gm3(m13+m23)(m1+m2)7/2(1e22)7/2a17/2a25e22 (33)
Ap6 = 47251024Gm3(m15+m25)(m1+m2)11/2(1e22)11/2a111/2a27e22(2+e22) (34)
Ap8 = 11025262144Gm3(m17+m27)(m1+m2)15/2(1e22)15/2a115/2a29e22(528+880e22+165e24) (35)
Ap10 = 2182958388608Gm3(m19+m29)(m1+m2)19/2(1e22)19/2a119/2a211e22(1664+5824e22+3640e24+364e26) (36)
B=Bp2+Bp4+Bp6+Bp8+Bp10 (37)
Bp2 = 34Gm3(m1+m2)1/2(1e22)3/2a13/2a23 (38)
Bp4 = 45128Gm3(m13+m23)(m1+m2)7/2(1e22)7/2a17/2a25(4+13e22) (39)
Bp6 = 5252048Gm3(m15+m25)(m1+m2)11/2(1e22)11/2a111/2a27(8+76e22+33e24) (40)
Bp8 = 11025262144Gm3(m17+m27)(m1+m2)15/2(1e22)15/2a115/2a29(64+1200e22+1720e24+305e26) (41)
Bp10 = 2182958388608Gm3(m19+m29)(m1+m2)19/2(1e22)19/2a119/2a211(128+3968e22+11872e24+7000e26+679e28) (42)
C=Cp3+Cp5+Cp7+Cp9+Cp11 (44)
Cp3 = 1516Gm3(m1m2)(m1+m2)3/2(1e22)5/2a15/2a24e2 (45)
Cp5 = 105256Gm3(m14m24)(m1+m2)9/2(1e22)9/2a19/2a26(4+3e22)e2 (46)
Cp7 = 472516384Gm3(m16m26)(m1+m2)13/2(1e22)13/2a113/2a28(8+20e22+5e24)e2 (47)
Cp9 = 24255524288Gm3(m18m28)(m1+m2)17/2(1e22)17/2a117/2a210(64+336e22+280e24+35e26)e2 (48)
Cp11 = 94594533554432Gm3(m110m210)(m1+m2)21/2(1e22)21/2a121/2a212(128+1152e22+2016e24+840e26+63e28)e2. (49)

حل المعادلتين (30) و(31) هو

x1 = C1cosB22ABt+C2sinB22ABtCBcosϖ2 (50)
y1 = 1B2Acos2ϖ2[(C2B22ABAC1sin2ϖ2)cosB22ABt(C1B22AB+ (51)
+AC2sin2ϖ2)sinB22ABt]CBsinϖ2,

حيث إن C1 وC2 ثابتان للتكامل. وبما أن مركبتي متجه الاختلاف المركزي صفريان في البداية، أي x10=0 وy10=0، فإن الثابتين C1 وC2 يُعطيان بـ

C1 = CBcosϖ2 (52)
C2 = CB(B2A)sinϖ2. (53)

وإذا عوضنا القيم أعلاه لثوابتنا في المعادلتين (50) و(51)، فإن الأخيرتين تصبحان:

x1 = CBcosϖ2cosB22ABt+CB22ABsinϖ2sinB22ABtCBcosϖ2 (54)
y1 = CBsinϖ2cosB22ABtCB22ABcosϖ2sinB22ABtCBsinϖ2. (55)

ومن ثم،

e1 = x12+y12 (56)
ϖ1 = arctan(y1x1). (57)

2.3 الحالة الخارجية

على غرار الحالة الداخلية، يمكننا اشتقاق الهاملتوني طويل الأمد للحالة التي يكون فيها الكوكب العملاق داخل مدار الكوكب الأصغر. ويُحصل على الهاملتوني في هذه الحالة بطريقة مشابهة للحالة الداخلية؛ والفرق الوحيد هو أننا نهمل الآن حدود O(e23) وما فوقها بدلا من حدود لا تقل رتبتها عن O(e13). ونتيجة لذلك، يكون الهاملتوني الاضطرابي طويل الأمد:

Hp2 = 14𝒢m0m1m2(m0+m1)(1e22)3/2a12a23(1+32e12) (58)
Hp3 = 1516𝒢m0m1m2(m0m1)(m0+m1)2(1e22)5/2a13a24e1e2(1+34e12)cos(ϖ1ϖ2) (59)
Hp4 = 964𝒢m0m1m2(m03+m13)(m0+m1)4(1e22)7/2a14a25[1+5e12+158e14+(32+152e12+4516e14)e22+(354+358e12)e12e22× (60)
×cos(2ϖ12ϖ2)]
Hp5 = 10564𝒢m0m1m2(m04m14)(m0+m1)5(1e22)9/2a15a26(1+52e12+58e14)e1e2cos(ϖ1ϖ2) (61)
Hp6 = 25256𝒢m0m1m2(m05+m15)(m0+m1)6(1e22)11/2a16a27[1+212e12+1058e14+3516e16+(5+1052e12+5258e14+ (62)
+17516e16)e22+(1894+3154e12+94564e14)e12e22cos(2ϖ12ϖ2)]
Hp7 = 47252048𝒢m0m1m2(m06m16)(m0+m1)7(1e22)13/2a17a28(1+214e12+358e14+3564e16)e1e2cos(ϖ1ϖ2) (63)
Hp8 = 12258192𝒢m0m1m2(m07+m17)(m0+m1)8(1e22)15/2a18a29[12+9e12+1898e14+1058e16+315256e18+(214+ (64)
+1892e12+396916e14+220516e16+6615512e18)e22+(2974+20798e12+1039564e14+
+2079128e16)e12e22cos(2ϖ12ϖ2)]
Hp9 = 242558192𝒢m0m1m2(m08m18)(m0+m1)9(1e22)17/2a19a210(1+9e12+634e14+10516e16)e1e2cos(ϖ1ϖ2) (65)
Hp10 = 396932768𝒢m0m1m2(m09+m19)(m0+m1)10(1e22)19/2a110a211[12+554e12+4958e14+115516e16+5775256e18+(9+ (66)
+4952e12+44554e14+103958e16+51975128e18)e22+(7154+21452e12+4504532e14+
+1501532e16)e12e22cos(2ϖ12ϖ2)]
Hp11 = 945945262144𝒢m0m1m2(m010m110)(m0+m1)11(1e22)21/2a111a212(1+554e12+1654e14+115532e16)e1e2cos(ϖ1ϖ2). (67)

لاحظ أن Hp2 هو نفسه كما في الحالة الداخلية.

باتباع الإجراء نفسه كما سبق، ننتهي مرة أخرى إلى نظام من معادلتين تفاضليتين خطيتين غير متجانستين، لهما الشكل نفسه كما في الحالة الداخلية:

x˙2 = Asin2ϖ1x2+(B2Acos2ϖ1)y2+Csinϖ1 (68)
y˙2 = (B2Asin2ϖ1)x2Asin2ϖ1y2Ccosϖ1, (69)

مع

A=Ap4+Ap6+Ap8+Ap10 (70)
Ap4 = 315256GMm1m2(m13+m23)(m1+m2)5a14a211/2(2+e12)e12 (71)
Ap6 = 758192GMm1m2(m15+m25)(m1+m2)7a16a215/2(1008+1680e12+315e14)e12 (72)
Ap8 = 3675524288GMm1m2(m17+m27)(m1+m2)9a18a219/2(3168+11088e12+6930e14+693e16)e12 (73)
Ap10 = 39698388608GMm1m2(m19+m29)(m1+m2)11a110a223/2(91520+549120e12+720720e14+240240e16)e12 (74)
B=Bp2+Bp4+Bp6+Bp8+Bp10 (75)
Bp2 = 38GMm1m2(m1+m2)2a12a27/2(2+3e12) (76)
Bp4 = 45256GMm1m2(m13+m23)(m1+m2)5a14a211/2(8+54e12+22e14) (77)
Bp6 = 758192GMm1m2(m15+m25)(m1+m2)7a16a215/2(224+3360e12+4620e14+805e16) (78)
Bp8 = 3675524288GMm1m2(m17+m27)(m1+m2)9a18a219/2(384+10080e12+292332e14+17010e16+1638e18) (79)
Bp10 = 39698388608GMm1m2(m19+m29)(m1+m2)11a110a223/2(7040+285120e12+1420320e14+1737120e16+557865e18) (80)
C=Cp3+Cp5+Cp7+Cp9+Cp11 (81)
Cp3 = 1564GMm1m2(m1m2)(m1+m2)3a13a29/2(4+3e12)e1 (82)
Cp5 = 105512GMm1m2(m14m24)(m1+m2)6a15a213/2(8+20e12+5e14)e1 (83)
Cp7 = 4725131072GMm1m2(m16m26)(m1+m2)8a17a217/2(64+336e12+280e14+35e16)e1 (84)
Cp9 = 24255131072GMm1m2(m18m28)(m1+m2)10a19a221/2(16+144e12+252e14+105e16)e1 (85)
Cp11 = 9459458388608GMm1m2(m110m210)(m1+m2)12a111a225/2(32+440e12+1320e14+1155e16)e1. (86)

حل المعادلتين (68) و(69) له الشكل نفسه كما في الحالة الداخلية:

x2 = CBcosϖ1cosB22ABt+CB22ABsinϖ1sinB22ABtCBcosϖ1 (87)
y2 = CBsinϖ1cosB22ABtCB22ABcosϖ1sinB22ABtCBsinϖ1 (88)

و

e2 = x22+y22 (89)
ϖ2 = arctan(y2x2). (90)

وفي كلتا الحالتين الداخلية والخارجية، تكون فترة التذبذب

T=2πB22AB. (91)

2.4 الاختلاف المركزي الأعظمي والمتوسط

استنادا إلى النتائج التحليلية في الأقسام السابقة، يمكننا إيجاد تعبيرات للاختلاف المركزي الأعظمي ومتوسط مربع الاختلاف المركزي. وباستخدام المعادلات (54)، و(55)، و(87)، و(88) نحصل على:

ei2=xi2+yi2=(CB)2(1cosB22ABt)2+C2B(B2A)sin2B22ABt,i=1,2 (92)

كما يمكن ملاحظة ذلك، فإن النتيجة أعلاه مستقلة عن خط طول حضيض الكوكب الأكبر. لذلك، وبأخذ متوسط المعادلة أعلاه زمنيا، نحصل على:

ei2=32(CB)2+12C2B(B2A)=(CB)22B3AB2Ai=1,2. (93)

إضافة إلى ذلك، يمكن أن يبيَّن بسهولة من المعادلة (92) أن القيمة العظمى للاختلاف المركزي هي

eimax=2CBi=1,2. (94)

3 النتائج العددية

لاختبار صلاحية تقديراتنا التحليلية، أجرينا سلسلة من التجارب العددية تضمنت تكامل معادلات الحركة الكاملة ومعادلات الحركة الدهرية للهاملتوني الكامل حتى الرتبة الحادية عشرة. وبالنسبة إلى الأولى، استخدمنا المكامل السيمبلكتيكي ذي تحويل الزمن الذي طوّره Seppo Mikkola (Mikkola, 1997)، وهو مصمم خصيصا لتكامل الأنظمة الثلاثية الهرمية. ويستخدم الرمز إحداثيات جاكوبي القياسية، أي إنه يحسب متجهي الموضع والسرعة النسبيين للمدارين الداخلي والخارجي في كل خطوة زمنية. أما في التكامل العددي للمعادلات الدهرية، فقد استخدمنا مكامل Bulirsch-Stoer المعطى في Press et al. (1992).

في جميع المحاكاة كانت كتلة النجم مساوية لكتلة شمسية واحدة، بينما ضُبطت كتلة الكوكب الأصغر على كتلة أرضية واحدة. وبالنسبة إلى الكوكب الأكبر استخدمنا أربع كتل مختلفة: m2=10,1,0.1and0.03MJ (وهو يقارب 10M). وفي الحالة الداخلية ثُبت نصف المحور الرئيسي للكوكب الصغير عند 1 au، بينما تراوح نصف المحور الرئيسي للكوكب الأكبر بين 1.1 au و7 au، وأُعطي اختلافه المركزي قيما في المجال [0.05،0.4] بخطوة مقدارها 0.05. إن اختيار قيم أعلى لاختلاف مركزية الجسم المضطرب (>0.4) لن يكون ذا معنى في حالتنا، لأن التقديرات التحليلية الجديدة مهمة لنسب أصغر من أنصاف المحاور الرئيسية؛ فالنظام ذو اختلاف مركزية أعلى للكوكب العملاق سيتطلب نسبة أكبر لنصفي المحورين الرئيسيين كي يكون مستقرا. وفي هذه الحالة، يكون الحل القديم لـ Georgakarakos (2003) مع توسع كثيرات حدود لوجندر حتى P3 كافيا. وفي الدراسة الحالية، لكل زوج من (a2,e2)، أُعطيا خط طول الحضيض الابتدائي والشذوذ الحقيقي للكوكب الأكبر القيم 0,90,180,270، أي إننا حاكينا ستة عشر مدارا مختلفا. وأخيرا، كان الكوكب الأصغر يبدأ دائما على طول الاتجاه الموجب لمحور x. وبالمثل، في الحالة الخارجية، ثبّتنا نصف المحور الرئيسي للكوكب العملاق عند 0.5 au، بينما غُيّر نصف المحور الرئيسي للكوكب الأصغر من 0.7 au إلى 3 au. وضُبطت قيم اختلاف المركزية وخط طول الحضيض والشذوذ الحقيقي للكوكب العملاق والموضع الابتدائي للكوكب الأصغر بالطريقة نفسها كما في الحالة الداخلية. وفي جميع المحاكاة كان زمن التكامل مساويا لفترة دهرية تحليلية واحدة كما تُعطى في المعادلة (91).

كما ذكرنا في المقدمة، يهدف هذا العمل أساسا إلى أنظمة يكون فيها الكوكبان متقاربين نسبيا لكنهما لا يزالان في مدارات مستقرة (يمكن أخذ فكرة عن حد الاستقرار للأنظمة التي نبحثها من Georgakarakos, 2013; Petrovich, 2015)، لأن نماذجنا المطورة سابقا غير كافية لمثل هذه الأنظمة. وكانت الكميات التي اختبرناها هي متوسط مربع الاختلاف المركزي والاختلاف المركزي الأعظمي. ولكل من الحالتين الداخلية والخارجية وجدنا الخطأ النسبي بين المعادلة (93) ومتوسط مربع الاختلاف المركزي الناتج من التكامل العددي لمعادلات الحركة الكاملة. وفي الوقت نفسه، قسنا الكمية نفسها ولكن هذه المرة باستخدام معادلات الحركة الدهرية المشتقة من الهاملتوني الكامل (في الاختلافات المركزية) من الرتبة الحادية عشرة. وهذا يضمن أن الخطأ الصغير بين التقديرات التحليلية والحل العددي للمسألة الكاملة ليس نوعا من الأثر المصطنع. ولا يدل على أن التقديرات التحليلية تقريب جيد للمسألة الكاملة إلا خطأ صغير في الأولى مصحوب بخطأ صغير في الثانية. وطبقنا الإجراء نفسه على الاختلاف المركزي الأعظمي المحصل لزوج معين من (ai,ei)، (i=2,1 للحالتين الداخلية والخارجية على الترتيب). وأخيرا، لكل زوج (ai,ei)، سجلنا الاختلاف المركزي الأعظمي الذي بلغه الكوكب الأصغر أثناء التكامل العددي.

كانت التقديرات التحليلية عموما في توافق جيد مع نتائج المحاكاة العددية. وكما هو متوقع، تفشل التقديرات التحليلية عندما تصبح الاضطرابات قوية ويبلغ الاختلاف المركزي المتولد للكوكب الأصغر قيما متوسطة إلى عالية، أو عندما يقترب النظام من حد الاستقرار. ويحدث ذلك عادة عند قيم أعلى لاختلاف مركزية الكوكب العملاق ونسب أصغر لنصفي المحورين الرئيسيين. أما لبقية القيم في فضاء نصف المحور الرئيسي - الاختلاف المركزي، فقد تبين أن التقديرات التحليلية تقريب جيد للمسألة الكاملة. وكانت الأنظمة القريبة من رنينات الحركة المتوسطة استثناء من هذه القاعدة. ففي مثل هذه الحالة تفشل النظرية الدهرية، ولا يعمل النموذج النظري المطور هنا على نحو جيد. وتفشل الطريقة تماما قرب الموقع الدقيق لرنين الحركة المتوسطة، لكنها تتحسن تدريجيا كلما ابتعدنا عنه. وبطبيعة الحال، يعتمد مدى التأثير على كتلة الجسم المضطرب واختلاف مركزيته، إذ يعتمد عرض الرنين على هاتين المعلمتين (انظر مثلا Murray & Dermott, 1999).

تعرض الأشكال 2 و3 و4 بعض نتائجنا للحالتين الداخلية والخارجية لجسم مضطرب كتلته كتلة مشتري واحدة. ويبين الشكل 2 الأخطاء النسبية اللوغاريتمية في متوسط مربع الاختلاف المركزي والاختلاف المركزي الأعظمي، مصحوبة برسم يبين الاختلاف المركزي الأعظمي المتحقق. ويشير العمود الأيسر إلى الحالة الداخلية، بينما يشير العمود الأيمن إلى الحالة الخارجية. ويعرض الصف الأول الخطأ النسبي بين التقدير التحليلي ومتوسط مربع الاختلاف المركزي المستحصل من التكامل العددي لمعادلات الحركة الكاملة. وينبغي أن يؤخذ في الحسبان أن الخطأ النسبي في <ei2> قد يكون أكبر بعامل مقداره اثنان مقارنة بخطأ في متوسط الاختلاف المركزي. وفي الجانب الأيسر من الرسوم توجد مدارات غير مستقرة (هروب أو قذف الكوكب الأرضي) يُشار إليها بعدم وجود رمز، أو مدارات يكتسب فيها الاختلاف المركزي المضطرب الأعظمي قيما عالية. غير أننا كلما انتقلنا إلى يمين الرسوم أصبحت نسبة نصفي المحورين الرئيسيين أصغر وتحسنت نتائجنا التحليلية تحسنا ملحوظا.

عندما يصبح الاختلاف المركزي المضطرب الأعظمي أكبر من 0.2-0.25، تبدأ فرضية بقاء اختلاف مركزية الكوكب الأرضي صغيرا بالفشل، لأن حدود اختلاف مركزية الكوكب الأرضي من الرتبة الثالثة، المهملة في الهاملتوني الدهري، تصبح مهمة للتطور الديناميكي للنظام. لذلك يُتوقع أن تبدأ النظرية التحليلية بالفشل في مثل هذه الحالات. ومع ذلك، قد تُظهر النتائج التحليلية أحيانا توافقا جيدا مع النموذج الكامل لأنظمة يصبح فيها الاختلاف المركزي الأرضي أكبر من 0.25. ولسوء الحظ، قد يكون ذلك التوافق الجيد مصادفة، إذ إن إغفال حدود الرتبة الثالثة في الهاملتوني عُوّض بأثر آخر ما (مثل القرب من رنين حركة متوسطة). ولتجنب استخلاص نتيجة خاطئة بشأن دقة التقديرات التحليلية، أنشأنا رسوم الصف الثاني، التي تعرض الخطأ النسبي اللوغاريتمي في متوسط مربع الاختلاف المركزي بين التقديرات التحليلية ونتائج التكامل العددي لمعادلات الحركة الدهرية الكاملة من الرتبة الحادية عشرة. وعادة، عندما يصبح الاختلاف المركزي أكبر من 0.25، تؤدي حدود الاختلاف المركزي الأعلى رتبة، المحذوفة من الهاملتوني التقريبي، إلى حل تحليلي يبالغ في تقدير فترة العناصر المدارية الدهرية وسعتها. ولذلك لا يمكن قبول أن التقديرات النظرية تعمل جيدا إلا عندما يكون كل من الخطأ بين التقديرات التحليلية ونتائج المعادلات الكاملة، والخطأ بين التقديرات التحليلية ونتائج المعادلات الدهرية الكاملة، صغيرا. والفجوة في الزاوية العلوية اليسرى من تلك الرسوم سببها أن الاختلاف المركزي المضطرب نما إلى أكبر من واحد، فلم يعد بالإمكان تعريف الكمية (1ei2)، التي تظهر في بعض مقامات معادلات الحركة الدهرية، ضمن مجموعة الأعداد الحقيقية، ومن ثم توقفت المحاكاة.

يقارن الصف الثالث من الرسوم التقدير التحليلي للاختلاف المركزي الأعظمي بالاختلاف المركزي الأعظمي الذي تحقق أثناء التكامل العددي لمعادلات الحركة الكاملة. ومرة أخرى، يكون الخطأ كبيرا للأنظمة في المنطقة العلوية اليسرى من الرسوم، ثم يصغر تدريجيا كلما انتقلنا إلى اليمين. ويعرض الصف الرابع من الرسوم الخطأ اللوغاريتمي النسبي بين التقدير التحليلي للاختلاف المركزي الأعظمي والاختلاف المركزي الأعظمي الناتج أثناء التكامل العددي لمعادلات الحركة الدهرية الكاملة من الرتبة الحادية عشرة. وكانت أسباب إدراج هذا الصف هي نفسها كما في حالة متوسط مربع الاختلاف المركزي. ويعرض الصف الأخير الاختلاف المركزي العددي الأعظمي الذي حصلنا عليه عند تكامل معادلات الحركة الكاملة. كما أن استخدام هذه الرسوم لتحديد مناطق الاختلاف المركزي العالي مؤشر مفيد على المواضع التي قد تواجه فيها تقديراتنا التحليلية مشكلات.

كما ذكرنا سابقا، عندما تكون الشروط الابتدائية لنظام ما في جوار رنين حركة متوسطة، تفشل نظريتنا. وكلما ابتعدنا عن تلك المنطقة تتحسن النتائج من جديد. وفي رسوم الصف الأول في الشكل 2، أشرنا إلى بضعة مواقع لرنين الحركة المتوسطة تتطابق مع أخذ عينات نصف المحور الرئيسي الابتدائي. ويتجلى الأثر الرنيني بخط عمودي أكثر حمرة (دلالة على تغير مفاجئ في الخطأ) مقارنة بالمنطقة المجاورة. وكلما انتقلنا إلى نسب أصغر لنصفي المحورين الرئيسيين، أصبح تجنب الرنينات أكثر صعوبة لأنها تصبح أقوى وتغطي مساحات أكبر من فضاء المعلمات. وهذا، إلى جانب حقيقة أن تقارب توسع متسلسلة كثيرات حدود لوجندر يصبح أبطأ كثيرا عند نسب أصغر لنصفي المحورين الرئيسيين، يجعل الأمور عسيرة إلى حد ما على تقديراتنا عندما a2/a11.9.

يبين الشكلان 3 و4 تحسن التقديرات التحليلية للحالتين الداخلية والخارجية كلما أضفنا حدودا أكثر إلى التوسع. وهما يعرضان الخطأ النسبي المئوي في الفترة وفي الاختلاف المركزي الأعظمي بين حلول تحليلية مقدرة عند رتب مختلفة من نسبة نصفي المحورين الرئيسيين. وبصورة أدق، حسبنا الكمية 100(K3Ki)/Ki مع i=5,7,9,11، حيث إن K إما فترة تطور الاختلاف المركزي/خط طول الحضيض، المعطاة بالمعادلة (91)، أو الاختلاف المركزي الأعظمي، وتشير الأدلة إلى رتبة التوسع التي استخدمناها. وتكون الزيادة في رتبة التوسع بين حل تحليلي وتحسينه الأول دائما مقدارها اثنان، لأن حدي A وB في الحل التحليلي يأتيان من كثيرات حدود لوجندر ذات أدلة زوجية، بينما يحتوي الحد C فقط على كثيرات حدود لوجندر ذات رتب فردية. وينبغي للقارئ أن يضع في اعتباره أن التغيرات في الحد A تؤثر في الفترة فقط، وأن التغيرات في الحد B تؤثر في كل من الفترة والقيمة العظمى للاختلاف المركزي، وأن التغيرات في الحد C تؤثر في الاختلاف المركزي الأعظمي فقط. لذلك فمن الضروري أن نرتقي في التوسع برتبتين في كل مرة للحصول على الصورة الصحيحة.

ما نجده من فحص الشكلين 3 و4 هو أن التوسع عالي الرتبة مهم للحصول على صورة أدق لتطور النظام مع تناقص نسبة نصفي المحورين الرئيسيين a2/a1. ويكون أثر تضمين الحدود الأعلى رتبة في الحل التحليلي أوضح بكثير في الفترة منه في القيمة العظمى للاختلاف المركزي. إضافة إلى ذلك، يزداد أثر تضمين الحدود الأعلى رتبة مع زيادة اختلاف مركزية الجسم المضطرب. وفي جميع الحالات، يكون التوسع حتى الرتبة الثالثة غير كاف للوصف الديناميكي للنظام، ويكون استخدام توسع لا يقل عن الرتبة الخامسة ضروريا، حتى في الأنظمة ذات نسب أنصاف محاور رئيسية أكبر. وأخيرا، نود الإشارة إلى أن الرسوم في الشكلين 3 و4 لا تتأثر بتغير كتلة الجسم المضطرب. وهذا متوقع كما توحي المعادلات المشتقة في الأقسام السابقة.

يمكن العثور على مزيد من نتائج محاكاتنا في الملحق C، حيث يمكن للقارئ أن يجد رسوما مشابهة لما في الشكل 2، ولكن لجسم مضطرب كتلته 10 و0.1 و0.03 من كتل المشتري. وتظهر النتائج الخاصة بتلك الكتل الأخرى سلوكا مماثلا لنتائج نظام كان فيه للجسم المضطرب كتلة مقدارها 1MJ. والفرق الرئيس هو أنه كلما صغرت (كبرت) كتلة الجسم المضطرب، ومن ثم ضعفت (قويت) الاضطرابات، أصبحت تقديراتنا التحليلية تقريبا جيدا للمسألة الكاملة على مساحة أكبر (أصغر) من فضاء معلماتنا.

Refer to caption
Refer to caption
Figure 2: الأخطاء في متوسط مربع الاختلاف المركزي والاختلاف المركزي الأعظمي لنظام ذي m0=1M,m1=1M,m2=1MJ (العمود الأيسر-الحالة الداخلية) ونظام ذي m0=1M,m1=1MJ,m2=M (العمود الأيمن-الحالة الخارجية). لكلا العمودين، من الأعلى: يعرض الصف الأول الخطأ النسبي اللوغاريتمي في متوسط مربع الاختلاف المركزي بين نتائج التكامل العددي لمعادلات الحركة الكاملة وتقديراتنا التحليلية، أي ri=log10[(ein2ei2)/ein2]. ويعرض الصف الثاني الخطأ النسبي اللوغاريتمي في متوسط مربع الاختلاف المركزي بين نتائج التكامل العددي لمعادلات الحركة الدهرية الكاملة من الرتبة الحادية عشرة وتقديراتنا التحليلية، أي ris=log10[(eis2ei2)/eis2]. ويعرض الصف الثالث الخطأ النسبي اللوغاريتمي في الاختلاف المركزي الأعظمي بين نتائج التكامل العددي لمعادلات الحركة الكاملة وتقديراتنا التحليلية، أي rim=log10[(einmaxeimax)/einmax]. ويعرض الصف الرابع الخطأ النسبي اللوغاريتمي في الاختلاف المركزي الأعظمي بين نتائج التكامل العددي لمعادلات الحركة الدهرية الكاملة من الرتبة الحادية عشرة وتقديراتنا التحليلية، أي rims=log10[(eismaxeimax)/eismax]. وأخيرا، يعرض الصف الخامس الاختلاف المركزي الأعظمي للكوكب الأصغر einmax الذي تحقق أثناء التكامل العددي لمعادلات الحركة الكاملة. في الحالة الداخلية i=1، بينما في الحالة الخارجية i=2. ويشير الدليل n إلى التكامل العددي لمعادلات الحركة الكاملة.
Refer to caption
Refer to caption
Refer to caption
Refer to caption
Figure 3: الخطأ المئوي للفترة T والاختلاف المركزي الأعظمي emax بين حلول تحليلية مبنية على رتب توسع مختلفة لنظام ذي m0=1M,m1=1M,m2=1MJ. في رسوم العمود الأيسر لدينا e2=0.1، بينما في رسوم العمود الأيمن لدينا e2=0.3. ومنحنيات الخطأ المعروضة هنا ناتجة من المقارنة بين النماذج التحليلية الآتية: P3P5، P3P7، P3P9، P3P11، حيث يشير Pn (n=3,5,7,9,11) إلى حل تحليلي مبني على توسع لكثيرات حدود لوجندر حتى الرتبة n.
Refer to caption
Refer to caption
Refer to caption
Refer to caption
Figure 4: الخطأ المئوي للفترة T والاختلاف المركزي الأعظمي emax بين حلول تحليلية مبنية على رتب توسع مختلفة لنظام ذي m0=1M,m1=1MJ,m2=1M. في رسوم العمود الأيسر لدينا e1=0.1، بينما في رسوم العمود الأيمن لدينا e1=0.3. ومنحنيات الخطأ المعروضة هنا ناتجة من المقارنة بين النماذج التحليلية الآتية: P3P5، P3P7، P3P9، P3P11، حيث يشير Pn (n=3,5,7,9,11) إلى حل تحليلي مبني على توسع لكثيرات حدود لوجندر حتى الرتبة n.

4 الخلاصة والمناقشة

حصلنا في هذا العمل على حل تحليلي للحركة طويلة الأمد لكوكب أصغر تحت تأثير الاضطرابات الثقالية لكوكب أكبر، مع وقوع جميع الأجسام في المستوى نفسه. وكان الكوكب الأصغر في البداية على مدار دائري. غير أن النظرية المطورة في هذه الورقة يمكن أن تنطبق على أي نظام يكون للكوكب الأرضي فيه اختلاف مركزي ابتدائي صغير غير صفري. وتُدخل الشروط الابتدائية المختلفة في معادلاتنا عبر الثابتين C1 وC2.

أجرينا توسعا عالي الرتبة للهاملتوني الاضطرابي بدلالة نسبة نصفي المحورين الرئيسيين للكوكبين، وحللنا معادلات الحركة المشتقة تحليليا. وقد قُدم التوسع والحل التحليلي لمعادلات الحركة بطريقة تتيح، تبعا للتطبيق والدقة المطلوبة، الاحتفاظ بالرتبة المرغوبة في نسبة نصفي المحورين الرئيسيين واختلاف مركزية الجسم المضطرب. علاوة على ذلك، يمكن دمج الحل الدهري الجديد بسهولة، عند الحاجة، مع الحدود قصيرة الفترة (كما في المثال Georgakarakos, 2003). وتصبح مثل هذه الحدود ضرورية في أنواع الأنظمة التي نبحثها عندما يكون اختلاف مركزية الجسم المضطرب منخفضا (0.01). ثم قورن الحل التحليلي بنتائج مستحصلة من محاكاة عددية، وكان الاتفاق بين النتائج التحليلية والعددية جيدا عموما.

وكما هو متوقع، لا تعمل النتائج التحليلية الجديدة جيدا قرب رنينات الحركة المتوسطة، وعندما تكون كتلتَا الجسمين الكوكبيين متقاربتين، ومن ثم لا يمكن افتراض ثبات اختلاف مركزية الجسم المضطرب وحضيضه بعد ذلك، وفي الحالات التي يتغير فيها حضيض الجسم الأكبر مع الزمن بسبب النسبية العامة. وينبغي معالجة مثل هذه الحالات بطريقة مختلفة. فعلى سبيل المثال، في الحالة التي يكون فيها مدار الكوكب العملاق داخل مدار الكوكب الأصغر، قد يلزم تضمين تصحيح ما بعد نيوتني لحركة حضيض الكوكب العملاق من أجل نمذجة النظام بدقة أكبر. وفي مثل هذا السيناريو، لن تعود للمعادلتين (68) و(69) معاملات ثابتة، وينبغي البحث عن نوع مختلف من الحلول. ومن ناحية أخرى، على الرغم من أن الكواكب القريبة من رنين حركة متوسطة يمكن أن تفسد حلنا تماما، فإننا كلما ابتعدنا عن موقع رنين الحركة المتوسطة، قد يوفر الحل الدهري تقريبا معقولا لحركة النظام، لكنه قد يحتاج إلى الاقتران بنمذجة ما لرنين الحركة المتوسطة لتحسين دقته. وأخيرا، ستزداد المسألة تعقيدا إذا افترضنا أن خط طول حضيض الجسم المضطرب يتأثر بالجسم المضطَرَب. وهذا افتراض ضروري إلى حد ما عند التعامل مع كوكبين متقاربين في الحجم، وقد أكدت ذلك بعض محاكاة الاختبار التي أجريناها. وفي هذه الحالة، وعلى خلاف ما يحدث عند تضمين التصحيح ما بعد النيوتني الذي ينتج عادة معدل تغير ثابت لحركة الحضيض، فإن معدل تغير الحضيض سيكون دالة غير ثابتة في الزمن، وقد يكون إيجاد حل لتلك المعادلات التفاضلية صعبا.

كما ذُكر سابقا، ارتفع عدد الكواكب الخارجية المكتشفة مؤخرا إلى أكثر من ثلاثة آلاف. وقد تحسنت طرائق الكشف بما يكفي على مدى السنين بحيث أصبحنا قادرين حاليا على كشف كواكب ذات حجم مقارب لحجم الأرض، وإن كانت أقرب بكثير إلى النجم المضيف (مثلا Marcy et al., 2014). وقد أصبح المجتمع العلمي أكثر اهتماما من أي وقت مضى بالإجابة عن سؤال ما إذا كانت الحياة كما نعرفها يمكن أن تزدهر في مواضع أخرى من الكون، ونتيجة لذلك ظهر عدد كبير من الأوراق التي تناقش قابلية السكن. ويعد الاختلاف المركزي الكوكبي جانبا مهما من قابلية السكن، لأنه يمكن أن يؤثر في مناخ الكوكب (e.g. Williams & Pollard, 2002; Spiegel et al., 2010; Dressing et al., 2010; Linsenmeier et al., 2015; Eggl et al., 2012). ولذلك، يمكن استخدام التقديرات التحليلية الجديدة المعروضة في الأقسام السابقة في مثل هذه الدراسات. وتشمل خططنا المستقبلية تحسين الحل للحركة طويلة الأمد بأخذ تأثيرات ديناميكية أخرى في الاعتبار، مثل رنينات الحركة المتوسطة أو النسبية العامة، ومحاولة تقديم حلول مدارية لأنظمة ذات كتل كوكبية متقاربة.

الشكر والتقدير

نود أن نشكر المحكم المجهول على تعليقاته المفيدة التي ساعدتنا في تحسين المخطوطة. كما نود أن نشكر فريق موارد الحوسبة عالية الأداء في New York University Abu Dhabi، وبخاصة Jorge Naranjo، لمساعدتنا في محاكاتنا العددية.

References

  • Brouwer & Clemence (1961) Brouwer D., Clemence G. M., 1961, Methods of celestial mechanics
  • Dressing et al. (2010) Dressing C. D., Spiegel D. S., Scharf C. A., Menou K., Raymond S. N., 2010, ApJ, 721, 1295
  • Eggl et al. (2012) Eggl S., Pilat-Lohinger E., Georgakarakos N., Gyergyovits M., Funk B., 2012, ApJ, 752, 74
  • Ford et al. (2000) Ford E. B., Kozinsky B., Rasio F. A., 2000, ApJ, 535, 385
  • Georgakarakos (2002) Georgakarakos N., 2002, MNRAS, 337, 559
  • Georgakarakos (2003) Georgakarakos N., 2003, MNRAS, 345, 340
  • Georgakarakos (2004) Georgakarakos N., 2004, Celestial Mechanics and Dynamical Astronomy, 89, 63
  • Georgakarakos (2006) Georgakarakos N., 2006, MNRAS, 366, 566
  • Georgakarakos (2009) Georgakarakos N., 2009, MNRAS, 392, 1253
  • Georgakarakos (2013) Georgakarakos N., 2013, New Astron., 23, 41
  • Harrington (1968) Harrington R. S., 1968, AJ, 73, 190
  • Laskar & Boué (2010) Laskar J., Boué G., 2010, A&A, 522, A60
  • Lee & Peale (2003) Lee M. H., Peale S. J., 2003, ApJ, 592, 1201
  • Linsenmeier et al. (2015) Linsenmeier M., Pascale S., Lucarini V., 2015, Planet. Space Sci., 105, 43
  • Marcy et al. (2014) Marcy G. W., et al., 2014, ApJS, 210, 20
  • McArthur et al. (2010) McArthur B. E., Benedict G. F., Barnes R., Martioli E., Korzennik S., Nelan E., Butler R. P., 2010, ApJ, 715, 1203
  • Migaszewski & Goździewski (2008) Migaszewski C., Goździewski K., 2008, MNRAS, 388, 789
  • Mikkola (1997) Mikkola S., 1997, Celestial Mechanics and Dynamical Astronomy, 67, 145
  • Murray & Dermott (1999) Murray C. D., Dermott S. F., 1999, Solar system dynamics
  • Naoz et al. (2013) Naoz S., Farr W. M., Lithwick Y., Rasio F. A., Teyssandier J., 2013, MNRAS, 431, 2155
  • O’Toole et al. (2009) O’Toole S. J., Tinney C. G., Jones H. R. A., Butler R. P., Marcy G. W., Carter B., Bailey J., 2009, MNRAS, 392, 641
  • Ofir et al. (2014) Ofir A., Dreizler S., Zechmeister M., Husser T.-O., 2014, A&A, 561, A103
  • Ogihara et al. (2014) Ogihara M., Kobayashi H., Inutsuka S.-i., 2014, ApJ, 787, 172
  • Petrovich (2015) Petrovich C., 2015, ApJ, 808, 120
  • Press et al. (1992) Press W. H., Teukolsky S. A., Vetterling W. T., Flannery B. P., 1992, Numerical recipes in FORTRAN. The art of scientific computing
  • Raymond et al. (2006) Raymond S. N., Mandell A. M., Sigurdsson S., 2006, Science, 313, 1413
  • Spiegel et al. (2010) Spiegel D. S., Raymond S. N., Dressing C. D., Scharf C. A., Mitchell J. L., 2010, ApJ, 721, 1308
  • Williams & Pollard (2002) Williams D. M., Pollard D., 2002, International Journal of Astrobiology, 1, 61
  • Winn & Fabrycky (2015) Winn J. N., Fabrycky D. C., 2015, ARA&A, 53, 409

Appendix A معادلات الحركة

للتبسيط، أسقطنا جميع الأدلة.

بما أن G=L1e2، فإن لدينا:

dGdt = Hdg1e2dLdtLe1e2dedt=Hdgdedt=1e2Le(1e2dLdt+Hdg)
dedt = 1e2Le(1e2Hl+Hdg). (95)

وبما أن الهاملتوني مستقل عن l، إذن:

dedt=1e2LeHdg. (96)
dϖdt = HGdϖdt=HeeGdϖdt=1e2LeHe (97)

Appendix B الهاملتوني الاضطرابي الكامل من الرتبة الحادية عشرة

Hp2 = 14𝒢m0m1m2(m0+m1)(1e22)3/2a12a23(1+32e12) (98)
Hp3 = 1516𝒢m0m1m2(m0m1)(m0+m1)2(1e22)5/2a13a24e1e2(1+34e12)cos(ϖ1ϖ2) (99)
Hp4 = 964𝒢m0m1m2(m03+m13)(m0+m1)4(1e22)7/2a14a25[1+5e12+158e14+(32+152e12+4516e14)e22+(354+358e12)e12e22× (100)
×cos(2ϖ12ϖ2)]
Hp5 = 𝒢m0m1m2(m04m14)(m0+m1)5(1e22)9/2a15a26[[(315256e1+15752048e15+1575512e13)e23+(10564e1+525512e15+525128e13)e2]cos(ϖ1ϖ2)+ (101)
+(735512e13+22054096e15)e23cos(3ϖ13ϖ2)]
Hp6 = 𝒢m0m1m2(m05+m15)(m0+m1)6(1e22)11/2a16a27[[(2362532768e1678752048e1447252048e12)e24+(78751024e142362516384e1647251024e12)e22]× (102)
×cos(2ϖ12ϖ2)+(5197532768e143118565536e16)e24cos(4ϖ14ϖ2)+(37520483937516384e141312532768e1678754096e12)e24
26252048e148754096e16+(43754096e16125256131252048e142625512e12)e2225256525512e12]
Hp7 = 𝒢m0m1m2(m06m16)(m0+m1)7(1e22)13/2a17a28[[(2362516384e1+8268751048576e17+826875131072e15+49612565536e13)e25+(82687532768e15+49612516384e13+ (103)
+236254096e1+826875262144e17)e23+(165375131072e17+16537516384e15+992258192e13+47252048e1)e2]cos(ϖ1ϖ2)+[(14033252097152e17+
+2338875524288e15+467775131072e13)e25+(15592516384e13+77962565536e15+467775262144e17)e23]cos(3ϖ13ϖ2)+(8918912097152e17+
+891891524288e15)e25cos(5ϖ15ϖ2)]
Hp8 = 𝒢m0m1m2(m07+m17)(m0+m1)8(1e22)15/2a18a29[[(1273387516777216e18+636693758388608e16+127338751048576e14+1819125524288e12)e26+(42446251048576e18+ (104)
+60637532768e12+21223125524288e16+424462565536e14)e24+(25467751048576e18+254677565536e14+12733875524288e16+36382532768e12)e22]cos(2ϖ12ϖ2)+
+[(173423251048576e16+173423258388608e18+173423251048576e14)e24+(104053952097152e16+104053952097152e14+1040539516777216e18)e26]cos(4ϖ14ϖ2)+
+(150300158388608e16+644143516777216e18)e26cos(6ϖ16ϖ2)+(385875131072e12+1350562533554432e18+45018751048576e16+42875262144+81033751048576e14)e26+
+12862565536e16+3858752097152e18+110258192e12+(13505625524288e16+115762565536e12+4051687516777216e18+128625131072+24310125524288e14)e24+122516384+
+2572532768+81033754194304e18+4862025131072e14+2701125131072e16+23152516384e12)e22+23152565536e14]
Hp9 = 𝒢m0m1m2(m08m18)(m0+m1)9(1e22)17/2[[(5348227567108864e19+534822752097152e15+7640325524288e13+891371258388608e17+848925524288e1)e27+(53482275262144e15+ (105)
+891371251048576e17+84892565536e1+534822758388608e19+764032565536e13)e25+(458419532768e13+50935532768e1+32089365131072e15+320893654194304e19+
+53482275524288e17)e23+(15280651048576e19+152806532768e15+2182958192e13+2546775131072e17+242558192e1)e2]cos(ϖ1ϖ2)+[(2427925533554432e19+
+728377654194304e15+3468465524288e13+728377658388608e17)e27+(1213962752097152e17+404654258388608e19+1213962751048576e15+5780775131072e13)e25+
+(24279255524288e17+24279255262144e15+115615532768e13+80930852097152e19)e23]cos3ϖ13ϖ)+[(1932430533554432e19+450900458388608e17+
+270540274194304e15)e27+(450900452097152e17+193243058388608e19+270540271048576e15)e25]cos(5ϖ15ϖ2)+(156434858388608e17+
+46930455134217728e19)e27cos(7ϖ17ϖ2)]
Hp10 = 𝒢m0m1m2(m09+m19)(m0+m1)10(1e22)19/2a110a211[[(2979726751048576e14+993242252097152e12+208580872516777216e18+2085808725268435456e110+625742617516777216e16)e26+
+(417161745536870912e110+198648454194304e12+595945352097152e14+125148523533554432e16+41716174533554432e18)e28+(5959453516777216e110+851350565536e14+
+2837835131072e12+595945351048576e18+1787836051048576e16)e22+(41716174533554432e110+12514852352097152e16+59594535131072e14+19864845262144e12+ (106)
+4171617452097152e18)e24]cos(2ϖ12ϖ2)+[(1291214925134217728e18+18445927516777216e14+184459275268435456e110+77472895533554432e16)e28+
+(1844592752097152e14+18445927533554432e110+129121492516777216e18+7747289554194304e16)e24+(1844592752097152e14+18445927533554432e110+
+129121492516777216e18+7747289554194304e16)e26]cos(4ϖ14ϖ2)+[(19198822533554432e18+5759646751073741824e110+26878351533554432e16)e28+
+(44797252516777216e18+62716153516777216e16+1343917575536870912e110)e26]cos(6ϖ16ϖ2)+(13780488154294967296e110+
+41341464452147483648e18)e28cos(8ϖ18ϖ2)+(72201071251073741824e18+144402142567108864e16+12502358388608+61886632533554432e14+6876292516777216e12+
+8664128552147483648e110)e28+(22920975524288e12+28880428567108864e110+416745262144+2062887751048576e14+240670237533554432e18+4813404752097152e16)e26+
+396965536+218295131072e12+275051716777216e110+4584195524288e161964655262144e14+(17681895131072e14+247546538388608e110+41257755262144e16+196465565536e12+
+3572132768+2062887754194304e18)e22+(8664128552097152e16+51984771367108864e110+433206427533554432e18+3713197951048576e14+750141262144+
41257755524288e12)e24+229209758388608e18]
Hp11 = 𝒢m0m1m2(m010m110)(m0+m1)11(1e22)21/2a111a212[[(688316879254294967296e19+9833098275134217728e15+688316879251073741824e17+1376633758517179869184e111+3277699425134217728e13+ (107)
+5959453533554432e1)e29+(229438959752147483648e111+993242254194304e1+546283237516777216e13+1638849712516777216e15+114719479875536870912e19+
+114719479875134217728e17)e27+(98330982754194304e15+68831687925134217728e19+32776994254194304e13+6883168792533554432e17+13766337585536870912e111+
+595945351048576e1)e25+(8513505262144e1+4682427751048576e13+983309827533554432e19+98330982758388608e17+1966619655134217728e111+14047283251048576e15)e23+
+168719476736(8950304563200e17+247973806080e1+111878807040e111+10228919500800e15+2237576140800e19+
+3409639833600e13)e2]cos(ϖ1ϖ2)+[(64560746258589934592e111+71734162567108864e13+6456074625134217728e15+150641741251073741824e19+
+27115513425536870912e17)e29+(1936822387533554432e15+193682238752147483648e111+81346540275134217728e17+45192522375268435456e19+215202487516777216e13)e27+
+(1936822387516777216e15+193682238751073741824e111+45192522375134217728e19+8134654027567108864e1721520248758388608e13)e25+
+168719476736(6715957248000e13+30221807616000e15+8814693888000e19+31732897996800e17+
+472215744000e111)e23])cos(3ϖ13ϖ2)+[(15679038375536870912e172239862625134217728e15+1119931312517179869184e111+111993131251073741824e19)e29+
+(522634612533554432e15+36584422875134217728e17+26131730625268435456e19+261317306254294967296e111)e27+168719476736(22477469414400e17+
+8027667648000e19+501729228000e111+12844268236800e15)e25])cos(5ϖ15ϖ2)+[(10335366112517179869184e19+3445122037568719476736e111+
+206707322252147483648e17)e29+168719476736(3527804966400e17+2204878104000e19+183739842000e111)e27]cos(7ϖ17ϖ2)+
+168719476736(135763327700e19+20364499155e111)cos(9ϖ19ϖ2)e29]

Appendix C النتائج العددية

Refer to caption
Refer to caption
Figure 5: كما في الشكل 2، ولكن لنظام ذي m2=10MJ للحالة الداخلية (العمود الأيسر) وm1=10MJ للحالة الخارجية (العمود الأيمن).
Refer to caption
Refer to caption
Figure 6: كما في الشكل 2، ولكن لنظام ذي m2=0.1MJ للحالة الداخلية (الرسم الأيسر) وm1=0.1MJ للحالة الخارجية (الرسم الأيمن).
Refer to caption
Refer to caption
Figure 7: كما في الشكل 2، ولكن لنظام ذي m2=0.03MJ للحالة الداخلية (العمود الأيسر) وm1=0.03MJ للحالة الخارجية (العمود الأيمن).