Dynamical Stability of Planetary Systems

Examensarbete 30 hp May 2017

Dynamical Stability of Planetary Systems

Aikaterini Stergiopoulou

Masterprogrammet i fysik

Master Programme in Physics

Uppsala University

Department of Physics and Astronomy,

Astronomy and Space Physics.

Master Thesis:

Dynamical Stability of Planetary Systems

Aikaterini Stergiopoulou Supervisor: Ulrike Heiter

Uppsala, May 2017

Popul¨ arvetenskaplig Sammanfattning

I denna avhandling unders¨ oker vi den dynamiska stabiliteten hos planetsystem. Vi ger en kort beskrivning av banor och detekteringsmetoder f¨ or exoplaneter. D¨ arefter presenterar vi grunderna f¨ or stj¨ arnutveckling samt en matematisk beskrivning av den numeriska metoden som anv¨ ands f¨ or simulering av planetsystem.

Efter den teoretiska delen fokuserar vi f¨ orst p˚ a Kepler-11, som ¨ ar ett planetsystem med en stj¨ arna som liknar solen och sex planeter med banor n¨ ara varandra och n¨ ara stj¨ arnan. F¨ or den sj¨ atte planeten i detta system ¨ ar medelt¨ atheten, massan och excentriciteten inte k¨ anda. Vi k¨ or tre serier av simuleringar med olika kombinationer f¨ or dessa tre kvantiteter f¨ or att ta reda p˚ a hur systemet p˚ averkas n¨ ar vi ¨ andrar n˚ agra av deras ursprungliga f¨ orh˚ allanden. I 9 av 108 system

¨

overlever alla sex planeter fram till slutet av simuleringarna.

D¨ arefter implementerar vi f¨ or stj¨ arnan (Kepler-11) en konstant massf¨ orlustniv˚ a i 30 miljoner

˚ ar och k¨ or simuleringar f¨ or de nio stabila systemen fr˚ an tidigare k¨ orningarna. F¨ orutom Kepler- 11-systemet implementerar vi en konstant massf¨ orlusthastighet f¨ or solen och k¨ or simuleringar f¨ or nio konfigurationer av ett hypotetiskt system som endast best˚ ar av Solen, Jorden och Jupiter. I b˚ ada fallen, Kepler-11 och det hypotetiska solsystemet, flyttas planeterna ut˚ at n¨ ar stj¨ arnan f¨ orlorar massan. Att studera utvecklingen av planetbanor under de sena stadierna av stj¨ arnutvecklingen kan ge oss viktiga insikter om livsbetingelser i universum, eftersom den cirkumstell¨ ara beboeliga zonen ocks˚ a r¨ or sig ut˚ at. D¨ armed skulle en fr˚ an b¨ orjan beboelig planet kunna forts¨ atta att vara v¨ ard f¨ or liv ¨ aven n¨ ar stj¨ arnan blir en r¨ od j¨ atte.

1

Abstract

The study of dynamical stability in planetary systems has become possible during the last few decades due to the development of numerical methods for long-term integrations of N-body systems. Since the 90’s the number of exoplanet detections has been increased significantly, making the simulations of other real planetary systems besides the Solar System feasible. One of the exciting new-found worlds is the system Kepler-11. Six planets which are located very close to each other orbit a solar-type star. In this project we first investigate the behavior of Kepler-11 when we change some of the initial conditions of the outermost planet of the system and then we approximate the Red Giant phase of solar-type stars in order to see how the planetary orbits are altered. For the first part we run three series of simulations (groups A,B,C). Each group has a different value for the mean density of planet Kepler-11g (1.0,1.5,2.0 g/cm ³ ). We run simulations for 36 different combinations of mass and eccentricity of planet Kepler-11g for each group. In nine configurations all six planets of the system continue to orbit the star until the end of the simulations. These nine stable configurations of Kepler-11 are used in the second part where we implement a constant mass- loss rate for the star which results in 30% mass loss after 30 million years, trying to approximate that way the mass loss of solar-type stars in Red Giant Branch. We also run nine simulations of a hypothetical system consisting only of the Sun, Earth and Jupiter where we implement the constant mass-loss rate to the Sun. In the Kepler-11 system, the orbits of planets Kepler-11g and Kepler-11e change by ∼45% and ∼54% respectively, after 30 million years, due to the mass loss of the star, while in the hypothetical planetary system the orbits of the two planets change by ∼43%.

The study of orbits and how they move outward during the Post-Main Sequence evolution of stars is essential for our understanding of the existence of a Habitable Zone, not just around stars in Main-Sequence phase, but also around stars in late stages of their evolution.

2

in Memory of Kiriakos Stergiopoulos,

my beloved father

CONTENTS 4

Contents

1 Introduction 12

2 Orbits 13

2.1 Orbital Elements . . . . 14

2.2 Kepler’s third law . . . . 14

2.3 Orbital Resonances . . . . 14

3 Types of exoplanets 15 4 The Habitable Zone 17 5 Observational techniques for the detection of extrasolar planets 17 5.1 Radial Velocity . . . . 17

5.2 Astrometry . . . . 18

5.3 Timing . . . . 19

5.4 Microlensing . . . . 19

5.5 Direct Imaging . . . . 20

5.6 Transits . . . . 21

6 Evolution of solar-type stars 23 6.1 The H-R diagram . . . . 23

6.2 The Main-Sequence (MS) . . . . 25

6.3 The Red Giant Branch (RGB) . . . . 27

6.4 Mass-loss . . . . 28

7 Symplectic Integrators for N-Body Problems 29 7.1 Hamiltonian systems . . . . 29

7.2 Symplecticness . . . . 30

7.3 Symplectic Transformations . . . . 32

7.4 Symplectic Integrators . . . . 32

7.5 A Hybrid Symplectic N-Body Integrator (MERCURY)[9] . . . . 33

8 Kepler-11 34 8.1 Dynamical Stability of the system Kepler-11 for various configurations . . . . 34

8.1.1 Group A (d = 1.0g/cm ³ ) . . . . 35

8.1.2 Group B (d = 1.5g/cm ³ ) . . . . 38

8.1.3 Group C (d = 2.0g/cm ³ ) . . . . 42

8.2 Orbital evolution during the Red Giant Branch(RG phase) . . . . 43

9 Discussion 47

10 Conclusion 49

CONTENTS 5

Appendices 51

A Kepler-11 for d=1.0g/cm ³ (group A) 51

A.1 8 Systems with only one planet left . . . . 51

A.2 24 Systems with two planets left . . . . 52

A.3 4 Systems with three planets left . . . . 55

A.4 Inclination of planet b for the first 500 years . . . . 56

A.5 Inclination of planet c for the first 500 years . . . . 57

A.6 Inclination of planet d for the first 2000 years . . . . 58

A.7 Inclination of planet e for the first 2000 years . . . . 59

A.8 Inclination of planet f for the first 1500 years . . . . 60

A.9 Inclination of planet g for the first 2000 years . . . . 61

B Kepler-11 for d=1.5g/cm ³ (group B) 62 B.1 8 Systems with three planets left . . . . 62

B.2 13 Systems with four planets left . . . . 63

B.3 10 Systems with five planets left . . . . 65

B.4 Inclination of planet b for the first 15000 years . . . . 66

B.5 Inclination of planet c for the first 15000 years . . . . 67

B.6 Inclination of planet d for the first 40000 years . . . . 68

B.7 Inclination of planet e for the first 40000 years . . . . 69

B.8 Inclination of planet f for the first 20000 years . . . . 70

B.9 Inclination of planet g for the first 40000 years . . . . 71

B.10 Eccentricity of planet b for the first 800 years . . . . 72

B.11 Eccentricity of planet c for the first 800 years . . . . 73

B.12 Eccentricity of planet d for the first 2000 years . . . . 74

B.13 Eccentricity of planet e for the first 2000 years . . . . 75

B.14 Eccentricity of planet f for the first 2000 years . . . . 76

C Kepler-11 for d=2.0g/cm ³ (group C) 77 C.1 6 Systems with three planets left . . . . 77

C.2 9 Systems with four planets left . . . . 77

C.3 17 Systems with five planets left . . . . 78

D Mass-loss in Kepler-11 80 D.1 Group B . . . . 80

D.2 Group C . . . . 81

References 82

LIST OF FIGURES 6

List of Figures

1 An a-e diagram of the asteroids in the main belt. For the construction of the diagram,

the catalog of proper elements of 493574 asteroids from the AstDys-2 was used[4][19] . . . 15

2 Semi-major axis-mass diagram for 1324 exoplanets[5]. . . . 16

3 Semi-major axis-mass diagram for 1182 exoplanets [5]. . . . . 16

4 Semi-major axis-mass diagram for 1050 exoplanets (a < 6AU and m < 5M J )[5]. . . . 16

5 Semi-major axis-mass diagram for 252 exoplanets (a < 1.4AU and m < 0.1M _J )[5]. . . . . 16

6 Radial velocity measurements of HD 139357[12]. . . . 18

7 Radial velocity measurements of 42Dra[12]. . . . 18

8 Depiction of a transit of a planet in front of a star when it is viewed edge-on[44]. . . . 22

9 First, second, third and fourth contact of a transiting planet depicted by the numbers 1,2,3 and 4 respectively[39]. . . . 23

10 Transits of the six planets b,c,d,e,f and g of the system Kepler-11[23]. . . . 24

11 Hertzsprung-Russell diagram[1] . . . . 25

12 Evolutionary tracks of stars with different initial masses[3] . . . . 28

13 a of the planets left in eight systems . . . . 36

14 e of the planets left in eight systems . . . . 36

15 e of planets d, e and g for system 15 . . . . 38

16 i of planets d, e and g for system 15 . . . . 38

17 e of planets d, e and g for system 36 . . . . 38

18 i of planets d, e and g for system 36 . . . . 38

19 a-e diagram for system 15 . . . . 39

20 a-e diagram for system 22 . . . . 39

21 a-e diagram for system 34 . . . . 39

22 a-e diagram for system 36 . . . . 39

23 e of the six planets for system 11 . . . . 40

24 i of the six planets for system 11 . . . . 40

25 e of the six planets for system 27 . . . . 40

26 i of the six planets for system 27 . . . . 40

27 i of planet g for systems 1 and 36. . . . 42

28 a-e diagram of system 4 . . . . 42

29 a-e diagram of system 12 . . . . 42

30 e of the six planets for the system 11 . . . . 43

31 i of the six planets for the system 11 . . . . 43

32 a of system 11 of group B . . . . 45

33 a of system 33 of group B . . . . 45

34 a of system 11 of group C . . . . 45

35 a of system 33 of group C . . . . 45

36 a of Earth in system 1 . . . . 46

37 a of Earth in system 8 . . . . 46

38 a of Earth in systems 1 and 2 . . . . 47

39 a of Earth in system 8 . . . . 47

LIST OF FIGURES 7

40 a of Earth and Jupiter in system 9 . . . . 48

41 a of Jupiter in systems 1-9 . . . . 48

42 a of the terrestrial planets . . . . 48

43 a of Jupiter and Saturn . . . . 48

44 a of Uranus and Neptune . . . . 48

45 e of planets d and g of system 2 . . . . 54

46 i of planets d and g of system 2 . . . . 54

47 e of planets e and g of system 26 . . . . 54

48 i of planets e and g of system 26 . . . . 54

49 e of planets d and g of system 28 . . . . 54

50 i of planets d and g of system 28 . . . . 54

51 e of planets e and g of system 31 . . . . 55

52 i of planets e and g of system 31 . . . . 55

53 e of planets d and g of system 35 . . . . 55

54 i of planets d and g of system 35 . . . . 55

55 i of planet b . . . . 56

56 i of planet b . . . . 56

57 i of planet b . . . . 56

58 i of planet b . . . . 56

59 i of planet b . . . . 56

60 i of planet b . . . . 56

61 i of planet c . . . . 57

62 i of planet c . . . . 57

63 i of planet c . . . . 57

64 i of planet c . . . . 57

65 i of planet c . . . . 57

66 i of planet c . . . . 57

67 i of planet d . . . . 58

68 i of planet d . . . . 58

69 i of planet d . . . . 58

70 i of planet d . . . . 58

71 i of planet d . . . . 58

72 i of planet d . . . . 58

73 i of planet e . . . . 59

74 i of planet e . . . . 59

75 i of planet e . . . . 59

76 i of planet e . . . . 59

77 i of planet e . . . . 59

78 i of planet e . . . . 59

79 i of planet f . . . . 60

80 i of planet f . . . . 60

81 i of planet f . . . . 60

82 i of planet f . . . . 60

LIST OF FIGURES 8

83 i of planet f . . . . 60

84 i of planet f . . . . 60

85 i of planet g . . . . 61

86 i of planet g . . . . 61

87 i of planet g . . . . 61

88 i of planet g . . . . 61

89 i of planet g . . . . 61

90 i of planet g . . . . 61

91 e of planets c,e,g of system 3 . . . . 62

92 i of planets c,e,g of system 3 . . . . 62

93 e of planets d,e,g of system 21 . . . . 63

94 i of planets d,e,g of system 21 . . . . 63

95 e of planets c,d,e,g of system 5 . . . . 64

96 i of planets c,d,e,g of system 5 . . . . 64

97 e of planets c,d,e,g of system 16 . . . . 64

98 i of planets c,d,e,g of system 16 . . . . 64

99 e of planets b,d,e,f,g of system 26 . . . . 65

100 i of planets b,d,e,f,g of system 26 . . . . 65

101 i of planet b . . . . 66

102 i of planet b . . . . 66

103 i of planet b . . . . 66

104 i of planet b . . . . 66

105 i of planet b . . . . 66

106 i of planet b . . . . 66

107 i of planet c . . . . 67

108 i of planet c . . . . 67

109 i of planet c . . . . 67

110 i of planet c . . . . 67

111 i of planet c . . . . 67

112 i of planet c . . . . 67

113 i of planet d . . . . 68

114 i of planet d . . . . 68

115 i of planet d . . . . 68

116 i of planet d . . . . 68

117 i of planet d . . . . 68

118 i of planet d . . . . 68

119 i of planet e . . . . 69

120 i of planet e . . . . 69

121 i of planet e . . . . 69

122 i of planet e . . . . 69

123 i of planet e . . . . 69

124 i of planet e . . . . 69

125 i of planet f . . . . 70

LIST OF FIGURES 9

126 i of planet f . . . . 70

127 i of planet f . . . . 70

128 i of planet f . . . . 70

129 i of planet f . . . . 70

130 i of planet f . . . . 70

131 i of planet g . . . . 71

132 i of planet g . . . . 71

133 i of planet g . . . . 71

134 i of planet g . . . . 71

135 i of planet g . . . . 71

136 i of planet g . . . . 71

137 e of planet b . . . . 72

138 e of planet b . . . . 72

139 e of planet b . . . . 72

140 e of planet b . . . . 72

141 e of planet b . . . . 72

142 e of planet b . . . . 72

143 e of planet c . . . . 73

144 e of planet c . . . . 73

145 e of planet c . . . . 73

146 e of planet c . . . . 73

147 e of planet c . . . . 73

148 e of planet c . . . . 73

149 e of planet d . . . . 74

150 e of planet d . . . . 74

151 e of planet d . . . . 74

152 e of planet d . . . . 74

153 e of planet d . . . . 74

154 e of planet d . . . . 74

155 e of planet e . . . . 75

156 e of planet e . . . . 75

157 e of planet e . . . . 75

158 e of planet e . . . . 75

159 e of planet e . . . . 75

160 e of planet e . . . . 75

161 e of planet f . . . . 76

162 e of planet f . . . . 76

163 e of planet f . . . . 76

164 e of planet f . . . . 76

165 e of planet f . . . . 76

166 e of planet f . . . . 76

167 e of planets c,d,e,g of system 6 . . . . 78

168 i of planet c,d,e,g of system 6 . . . . 78

LIST OF FIGURES 10

169 e of planets c,d,e,f,g of system 22 . . . . 79

170 i of planet c,d,e,f,g of system 22 . . . . 79

171 e of planets c,d,e,f,g of system 35 . . . . 79

172 i of planet c,d,e,f,g of system 35 . . . . 79

173 a of system 6 from group B . . . . 80

174 a of system 15 from group B . . . . 80

175 a of system 27 from group B . . . . 81

176 a of system 31 from group C . . . . 81

177 a of system 34 from group C . . . . 81

LIST OF TABLES 11

List of Tables

1 Stellar properties of Kepler-11[17] . . . . 34

2 Planetary properties of the Kepler-11 system[17] . . . . 35

3 Group A simulations - d = 1.0g/cm ³ . . . . 37

4 Group B simulations - d = 1.5g/cm ³ . . . . 41

5 Group C simulations - d = 2.0g/cm ³ . . . . 44

6 Semi-major axes and eccentricities for Jupiter and Earth in nine configurations. . . . . . 46

7 Maximum values of Earth’s eccentricity and inclination for the first eight configurations. . 47

8 Final values of the semi-major axes of the planets of the Solar System. . . . . 49

9 Stable systems of groups B and C . . . . 50

10 Dynamical evolution of Kepler-11 for eight configurations (1,6,9,10,11,14,19,24 of Table 3) with only one planet left after the end of the simulations. Group A: d g = 1.0g/cm ³ . . . 51

11 Dynamical evolution of Kepler-11 for 24 configurations (2,3,4,5,7,8,12,13,16,17,18,20,21,23, 25,26,27,28,29,30,31,32,33,35 of Table 3) with two planets left after the end of the simu- lations. Group A: d g = 1.0g/cm ³ . . . . 52

13 Dynamical evolution of Kepler-11 for 4 conf. (15,22,34,36 of Table 3) with three planets left after the end of the simulations. Group A: d _g = 1.0g/cm ³ . . . . 55

14 Dynamical evolution of Kepler-11 for eight configurations (1,2,3,7,8,12,21,28 of Table 4) with three planets left after the end of the simulations. Group B: d _g = 1.5g/cm ³ . . . . . 62

15 Dynamical evolution of Kepler-11 for thirteen configurations (4,5,13,14,16,17,19,20,25,29, 31,35,36 of Table 4) with four planets left after the end of the simulations. Group B: d _g = 1.5g/cm ³ . . . . 63

17 Dynamical evolution of Kepler-11 for ten configurations (9,10,18,20,22,23,24,26,30,32,34 of Table 4) with five planets left after the end of the simulations. Group B: d g = 1.5g/cm ³ 65 18 Dynamical evolution of Kepler-11 for six configurations (1,7,8,14,21,24 of Table 5) with three planets left after the end of the simulations. Group C: d _g = 2.0g/cm ³ . . . . 77

19 Dynamical evolution of Kepler-11 for nine configurations (3,6,12,13,15,16,19,23,32 of Table 5) with four planets left after the end of the simulations. Group C: d _g = 2.0g/cm ³ . . . . 77

21 Dynamical evolution of Kepler-11 for 17 configurations (2,4,5,9,10,17,18,20,22,25,26,27,28,29, 30,35,36 of Table 5) with five planets left after the end of the simulations. Group C: d _g = 2.0g/cm ³ . . . . 78

23 Dynamical evolution of Kepler-11 for five configurations (6,11,15,27,33 of Table 4) of group B. . . . . 80

24 Dynamical evolution of Kepler-11 for four configurations (11,31,33,34 of Table 5) of group

C. . . . . 81

1 INTRODUCTION 12

1 Introduction

The first confirmed detection of exoplanets -planets of other planetary systems besides the Solar System- was made in 1992 by A. Wolszczan and D.A. Frail. They were searching for a stellar companion to the millisecond (6.2ms) radio pulsar PSR1257+12 using Arecibo radiotelescope, however they found two planets orbiting the pulsar and indications for the existence of a third planet[47].

Some other possible exoplanet detections, which were confirmed later, had been made even earlier[31].

In 1988, Campbell et al. suggested the existence of companions with masses in a range of 1 − 9M _J (where M _J is the mass of Jupiter) for seven stars, after measurements of their relative radial velocities[8]. In 1989, Latham et al. suggested, due to the detection of periodic variations in radial velocity, that an object of a mass of 11M _J orbits the star HD114762[22]. The first detection of a planet orbiting a solar- type star was made by Michel Mayor and Didier Queloz in 1995[28]. They measured the radial velocities of 142 stars and they found a hot Jupiter around 51 Pegasi[28]. After 1995 more planets were discovered by the radial velocity method[31]. The first planet detected by the gravitational microlensing method was confirmed in 2004[7][31]. The photometric transit technique was first used in 1999 for observations of already detected planets: G.W.Henry et al. and David Charbonneau et al. observed and measured the photometric transit of the planet HD 209458 b as it orbited the host star HD 209458. This planet had been already detected by the radial velocity method[15][10][31]. However, the first exoplanet discovery by the transit method was made in 2003 by Maciej Konacki et al.[20][31]. Today, there are about 3000 confirmed planets among which there are ∼ 500 multi-planet systems[17][5][2].

The goals of this project were to study and compare the different outcomes of the evolution of a planetary system when we change some of its initial conditions and to investigate the effects of mass loss during the RGB (Red Giant Branch) phase of a star to planetary orbits. To simulate the evolution of planetary systems, a hybrid symplectic N-body integrator was used[9]. The systems studied were Kepler-11 as well as several variations of the Solar System.

We perform 108 simulations which we divide in to 3 groups, A,B and C with a different value for the density of Kepler-11g, namely 1, 1.5 and 2.0g/cm ³ , respectively. Each group consists of 36 configurations with different combinations for the mass and the eccentricity of planet g. None of the systems in group A survive 8 million years intact. However, in 5 systems of group B and in 4 of group C all six planets continue orbiting the star until the end of the simulations(8Myr).

For the next part of the project we implement a constant mass-loss rate to the star which results in 30% mass loss after 30 million years, to approximate the behavior of solar-type stars while being in the RG phase. We investigate the evolution of the orbits in the stable configurations of Kepler-11 from groups B and C. We also study how mass-loss affects a hypothetical system consisting of only the Sun, Earth and Jupiter in which Jupiter’s orbit is located closer to the Earth’s one[43].

In Chapters 2,3 and 4 the basics about orbits, types of exoplanets and the Habitable Zone are

explained. In Chapter 5, the exoplanet detection methods are described. In Chapter 6, two stages of

stellar evolution namely the Main Sequence and the RG phase are analyzed, while in Chapter 7, the

theory of symplectic N-body integrators is presented. In Chapter 8 the results of the simulations are

shown. Then, a Discussion follows in Chapter 9 and finally, Chapter 10 concludes the report.

2 ORBITS 13

2 Orbits

For the following brief description of planetary orbits, the structure of chapter 2.1 of ”The Exoplanet Handbook” by Michael Perryman was used[31].

The geometrical shape of the path of a planet that orbits a star is an ellipse, as the first law of Kepler states. An ellipse, the center of which is located at the origin, can be described in Cartesian coordinates:

x ² a ² + y ²

b ² = 1 (1)

where a is the semi-major axis and b is the semi-minor axis of the ellipse, or in polar coordinates with the origin located in one of the foci:

r = a(1 − e ² )

1 + e cos θ (2)

where e is the eccentricity of the ellipse and it is related to a and b:

b ² = a ² (1 − e ² ) (3)

In the case of polar coordinates, where the star is located in one of the foci, the angular coordinate θ is called true anomaly. In order to find the pericentre (q) and the apocentre (Q) distance we substitute θ = 0 and θ = π respectively in eq.(2):

q = a(1 − e) (4)

Q = a(1 + e) (5)

Additionally to true anomaly, two more angles are used to describe the position of a celestial body in a particular moment. These are the eccentric (E) and the mean anomaly (M ). The relation between the eccentric and the true anomaly is[31]:

cos θ(t) = cos E(t) − e

1 − e cos E(t) (6)

The mean anomaly is defined as:

M = 2πt

P (7)

where P is the orbital period of the celestial body and t is the time measured after the object passes the pericentre.

Even though the angular rate of the motion of an orbiting object is not constant, an average angular rate can be defined and it is called mean motion:

n = 2π

P (8)

If we take into consideration the formula of mean motion the mean anomaly can be written:

M = nt (9)

2.1 Orbital Elements 14

The relation between the mean and the eccentric anomaly is known as Kepler’s relation[31]:

M (t) = E(t) − e sin E(t) (10)

2.1 Orbital Elements

In order to describe an orbit in three dimensions and to define where the object is at a particular moment, we need to specify 6 quantities which are called orbital elements. These are the semi-major axis a, eccentricity e, inclination i, argument of pericentre ω, longitude of the ascending node Ω and the mean anomaly M.

The semi-major axis and the eccentricity describe the elliptical orbit in two dimensions while the rest of the elements show us the orientation of the orbit and the current position of the object. However, sometimes the angles which are used differ. For example, in some cases is practical to use the longi- tude of pericentre p, the true longitude l and the mean longitude λ, which are related to the previous elements[31]:

p = Ω + ω (11)

l = p + θ (12)

λ = p + M (13)

2.2 Kepler’s third law

In the two-body problem, we can use Kepler’s third law to calculate the orbital period of an object when we know its semi-major axis and vice versa[31]. Kepler’s third law can be formulated either by considering the orbits of the two objects with respect to their common barycentre or by considering the orbit of the secondary object with respect to the central object. In the latter case, which is most often stated in textbooks, the form of the third law is:

P ² = 4πa ³

G(M S + M P ) (14)

where M S is the mass of the star or the central object, M P is the mass of the planet or the secondary object and a is the semi-major axis of the planet with respect to the star.

2.3 Orbital Resonances

The phenomenon of resonance has been observed and studied mainly in the Solar system[29]. A mean-

motion orbital resonance (MMR) occurs when the periods of two (or more) objects of a planetary

system, for example of two planets, a planet and a moon or a planet and an asteroid, follow a ratio

of small integers. Orbital resonances may either stabilize or disrupt a system. For example, the ratio

of Neptune’s and Pluto’s mean motions is 3:2. This orbit configuration prevents the two bodies from

3 TYPES OF EXOPLANETS 15

approaching each other, resulting in long-term stability[29]. Another example of resonances which seem necessary for the stability of the system occur between three of Jupiter’s moons, namely Io, Europa and Ganymede. The orbits of Io and Europa follow a 2:1 ratio while at the same time Europa and Ganymede follow also a 2:1 ratio. The ratio of their orbits is 4:2:1, which means that when Ganymede completes a revolution around Jupiter, Europa has completed two and Io four. This type of resonance is called Laplace resonance and prevents the conjunction of the three moons[29].

Examples of disruptive resonances can be found in the main asteroid belt between the orbits of asteroids and Jupiter. There are areas inside the main asteroid belt which are almost empty of asteroids, because Jupiter creates a dynamically-unstable environment at these spots due to orbital resonances.

The aforementioned areas are called Kirkwood gaps after the person who first observed them. In Fig.1 some of the Kirkwood gaps for specific values of the semi-major axis are visible: 2.06AU, 2.5AU, 2.82AU, 2.95AU, 3.27AU. These gaps correspond to 4:1, 3:1, 5:2, 7:3 and 2:1 mean motion ratios respectively.

0 0.1 0.2 0.3 0.4 0.5 0.6

1.5 2 2.5 3 3.5 4

e

a (AU)

Figure 1: An a-e diagram of the asteroids in the main belt. For the construction of the diagram, the catalog of proper elements of 493574 asteroids from the AstDys-2 was used[4][19]

3 Types of exoplanets

Among the confirmed extrasolar planets there are terrestrial rocky planets as well as giant gas ones that

resemble Jupiter (Jovian planets) and their mass is usually > 10M _E [31]. Planets with mass between 3

and 10 M _E are called either super earths or gas dwarfs when their composition does not consist mainly

3 TYPES OF EXOPLANETS 16

Figure 2: Semi-major axis-mass diagram for 1324 exoplanets[5].

Figure 3: Semi-major axis-mass diagram for 1182 exoplanets [5].

of solid material[31]. A subcategory of Jovian planets, is Hot Jupiters. They are gas giants that are located very close to the host star. Due to their short orbital period (P < 5 − 10d) and their large mass (m > 0.1M _J ) they are easily detected by radial velocity measurements and the probability of observing them transiting their host star is greater than for other planets[14][42].

In Fig.2,3,4 and 5, diagrams from the Extrasolar Planets Encyclopedia (exoplanet.eu[5]) of the semi- major axis versus the mass for a number of detected exoplanets are presented. In Fig.2 we see that the semi-major axes of the planets vary significantly. The majority of them seem to have semi-major axes either around the value 0.04-0.05AU or around 1-2AU. However, there are planets, with a up to 5000AU.

In Fig.3 we see the same diagram in a smaller scale, for a < 6AU and for m < 15M _J . A concentration of Jupiter-like planets (1 to 5M J ), the orbits of which are located very close to their star can be observed.

In Fig.4, in an even smaller scale (m < 5M _J ), hot Jupiters are easily noticeable. Finally, Fig.5 shows that super Earths close to their host star have been detected as well.

Figure 4: Semi-major axis-mass diagram for 1050 exoplanets

(a < 6AU and m < 5M J )[5].

Figure 5: Semi-major axis-mass diagram for 252 exoplanets

(a < 1.4AU and m < 0.1M _J )[5].

4 THE HABITABLE ZONE 17

4 The Habitable Zone

The human desire to search and find new worlds, new planets that resemble Earth -or not- and might be able to harbor life, is one of the major reasons that drives astronomers to look for extrasolar planets.

In order for carbon-based life to exist, the presence of liquid water is necessary.

Thus, the habitable zone (HZ) was defined as the region around a star where the planetary orbits are such, that water is sustained in liquid form on the surface of a planet for billions of years[31]. Some of the main factors that affect the location of the habitable zone in a planetary system are the spectral type of the star, the distance of the planet from the star, the eccentricity of the planet’s orbit, its rotation, its atmospheric properties, internal heating sources such as the existence of radionuclides, tidal heating and of course stellar evolution, either due to the increase of stellar luminosity during the Main Sequence phase or due to the change of stellar mass and radius when the star leaves the Main Sequence[18][31].

The inner and outer limits of the habitable zone of a star of 1M  in Main Sequence are 0.95AU and 1.37-1.67AU respectively according to a model proposed by Kasting et al. 1993[18]. Another, more recent model by Kopparapu et al.2013[21] suggests that the inner and outer limits are 0.99 and 1.70 AU respectively.

The habitable zone, however, is not stable through time. As the stellar luminosity changes and the star evolves based on its mass, the HZ moves outward[31]. The continuously habitable zone, that is the region, the limits of which remain inside the HZ at each given time, for the Solar System for the past 4Gyr are narrower (0.95-1.15AU) than the current HZ, as one might expect[31].

5 Observational techniques for the detection of extrasolar plan- ets

5.1 Radial Velocity

Measurements of the radial velocity of stars are used as an observational technique to detect exoplanets with which ∼620 planets have been confirmed so far[17].

The presence of a planet orbiting a star affects the motion of it. Both the planet and the star orbit their common barycenter. This revolution of the star around the barycenter causes changes to the expected spectra. As the star travels around the barycenter, the Doppler shift of the spectral lines is measured.

If r _c (t) is the distance of the star from the barycenter and assuming a coordinate system where the z-axis lies on the line of sight, the z-coordinate and the velocity projected onto the line of sight are given by[31]:

z(t) = r c (t) sin i sin (ω + θ) (15)

v r = sin i[ ˙r c (t) sin (ω + θ) + r(t) ˙ θ cos (ω + θ)] (16)

v _r = K[cos(ω + θ) + e cos ω] (17)

5.2 Astrometry 18

Figure 6: Radial velocity measurements of HD 139357[12].

Figure 7: Radial velocity measurements of 42Dra[12].

where K is the radial velocity semi-amplitude ((v _rmax − v _rmin )/2[26]) and a _c is the semi-major axis of the star with respect to the barycenter of the system:

K = 2π P

a _c sin i

(1 − e ² ) ^1/2 (18)

When the stellar mass is known, the measurements of radial velocity provide a lower limit for the mass of a planet M _lim = M _p sin i, which is determined by the factor sin i. In Fig. 6 and 7, two examples of radial velocity vs time (in days) measurements are presented. The points are the data and the solid line is the orbital solution[12].

5.2 Astrometry

One other observational technique which is used for the detection of planets is Astrometry. This method is based on the tracking and measurement of the precise positions and motions of stars which can point out the presence of an orbiting planet[31].

In the case of two bodies (a star and a planet), the semi-major axis of the star orbiting the center of mass of the system is given by:

a _c = M _P M S + M P

a ' M _P M S

a (19)

where a is the semi-major axis of the planet with respect to the star. a _c is a measurable quantity, thus it is preferable to be expressed in seconds of arc. The more planets a star hosts the more complex the motion of it with respect to the barycenter. The curves which depict the stellar motion around the barycenter, of a star that hosts a multi-planet system, are called planet mandalas after the Sanskrit word for circle[30][45].

In order to define the astrometric position of a star on the sky in terms of angular position, five

quantities are needed to be known: right ascension (a ₀ ), declination (δ ₀ ), the orthogonal components of

proper motion µ a? = µ a cos δ, µ δ and the parallax $[31].

5.3 Timing 19

One advantage of astrometry in comparison to radial velocity measurements is that instead of M _p sin i, the planetary mass M _p is defined due to determination of the inclination of the orbit through the astrometric measurements. However, in multiple planetary systems only the relative inclinations between two orbits can be defined[31]. So far, just one planet has been discovered by this method[17].

5.3 Timing

This method is based on the property of some types of primary objects to send out signals periodically, thus the existence of orbiting planets intercepts these signals and they can be detected. There are three types of primary objects which are used due to the periodicity of their signals for the detection of exoplanets, namely radio pulsars, pulsating stars and eclipsing binaries[31].

The first confirmed exoplanet discovery was made, as has already been mentioned, around the mil- lisecond radio pulsar PSR1257+12 where two terrestrial planets were found orbiting the pulsar[47]. The accuracy of timing pulsar radio emission that reaches our detectors allows the detection of planets that orbit pulsars. The timing of pulsars with spin periods of the order of 1 second can show the existence of Jovian and terrestrial planets, whereas the timing residuals of millisecond pulsars can indicate that even smaller objects, for example of the size of the moon or of a large asteroid, orbit the pulsar[31][46].

The amplitude of timing residuals in the case of a circular edge-on orbit of period P (in years) and a pulsar of mass 1.35M  is[46]:

τ p = 1.2 M _P M _Earth

!

P ^2/3 ms (20)

where M p and M Earth is the masses of the planet orbiting the pulsar and of Earth respectively. Five planets have been detected with pulsar timing variations[17]. The second type of primary objects that send out periodical signals are pulsating stars which are stars in their late stages of evolution. Only two planets have been detected so far with this method[17]. Finally, the existence of a planet around an eclipsing binary system which orbits both of the stars (P-type or circumbinary planets)[33] causes a deviation from the expected time of the primary and secondary eclipses between the two stars, because the binary stellar system follows an orbit around the common barycenter of the planetary system[35].

For timing accuracies of the order of 10 seconds, planets of mass of 10M _J with orbital periods of 10 to 20 years could be detected[35]. Eight planets have been detected with eclipse timing variation, with their orbital periods ranging from 416 days to ∼ 28 years[17].

5.4 Microlensing

According to general theory of relativity, light is deflected by mass distributions due to spacetime

distortion. Compact massive objects can bend the light in such a way that multiple rays reach the

observer from a single source making the source appear either as if it was multiple sources or resembling

a ring-shaped structure. When the source, the object that acts as lens and the observer are aligned and

the lens is axially symmetric then a ring-shaped image appears which is called Einstein ring. Otherwise,

when symmetry is broken, a multiple image is created (an image where more than one depictions of the

source appear)[6]. The angle of deflection of the light rays depends on the gravitational field of the mass

5.5 Direct Imaging 20

distribution[6]. This phenomenon allows us to observe objects which are located physically behind a massive object in our line of sight.

There are two types of gravitational lensing acknowledged so far, namely strong and weak[31]. When the ramifications of the lensing allow us to attribute intrinsic properties to individual objects which are part of the lensing, then the gravitational lensing is strong and results in multiple images and arcs[38].

Additionally, strong lensing can be further divided into macrolensing and microlensing depending on the resolution of the telescope[31]. On the other hand, in weak lensing the distortion and magnification of the source is small and it is not possible for individual sources to be traced, but an analysis can be performed in a statistical sense[38].

Gravitational microlensing can be considered as a case of strong lensing and it occurs when the separation of the expected images is too small to be resolved[31][38]. Then, the apparent brightness of the source changes according to the relative movement between the lens and the line of sight of the observer and the source[27][41]. The time scale of the duration of a microlensing event is between hours and years and some properties of this phenomenon are that the lightcurves are symmetric with respect to the maximum amplification of the brightness and that the events do not repeat[41]. Microlensing occurs for small compact masses in the range: 10 ⁻⁶ ≤ _M ^m

≤ 10 ⁶ [38].

When the object that acts as lens is a planetary system, consisting of a star and a planet (binary lens), the shape of the lightcurve differs from the case where a single star acts as lens and thus the existence of planets around stars can be detected[13]. There are three parameters one has to take into consideration in cases where the lens is binary, namely, the mass ratio ^M _M

S

, where M _P is the mass of the planet and M _S is the mass of the star, the projected distance between the planet and the star and the angle between the star-planet line and the relative to the lens path of the source[38]. The lightcurves of the source are affected the most when the projected distance s between the star and the planet is 0.6 ≤ _R ^s

E

≤ 1.6, where R _E is the Einstein radius, and this distance range is called lensing zone[38].

The lensing zone overlaps with the habitable zone making it possible for planets with liquid water to be detected with the microlensing method[38]. There have been ∼44 planets so far detected by this method[17].

5.5 Direct Imaging

This method refers to detecting planets as point sources either due to the reflected light from the host star or due to infrared radiation that the planet itself emits[31]. However, the main obstacle that one must overcome in order to be able to detect a planet is the brightness of the star. The ratio of planet to star brightness and their angular separation are two parameters of interest for the detection of exoplanets by imaging. The former depends on several properties of the system such as the stellar spectral type, the semi-major axis of the orbit of the planet, the planet’s mass and radius, the observation wavelength etc. The ratio of the planet to stellar flux for the stellar light that the planet reflects for observation wavelength λ can be written [31]:

f _p

f = p R _p a

! 2

g (21)

where p(λ) is the geometric albedo and g(α) is a phase dependent function. The values for extrasolar

5.6 Transits 21

planets are expected to be of the order of ∼ 10 ⁻¹⁰ for the reflected light (of solar type stars) and of

∼ 10 ⁻⁵ for the thermal emission[31]. The observed angular separation of a planet from its host star depends on the geometry of the planetary orbit and the distance of the system from us[31]. ∼44 planets have been confirmed by imaging and have orbital semi-major axes of tens to hundredths of AU[17].

5.6 Transits

The technique of finding transiting planets requires a specific geometry for the planetary systems with respect to our line of sight. When we observe a star and a planet passes in front of it, the planet obstructs a small part of the star decreasing this way the apparent brightness of it. Thus, it is impossible to detect exoplanets in ’face on’ planetary systems. The decrease of the star’s brightness is observed as a drop in its lightcurve which depends on the sizes of the star and the planet (Fig. 8). As of 2 Mar 2017, 2714 planets have been detected by this method making it the most effective technique so far[17].

The strategy which is followed for detecting transiting planets is the observation of large numbers of stars in wide-angle searches. The aim is the detection of the aforementioned decrease of the star’s brightness which repeats periodically[31]. This is a first indication for the existence of planets orbiting a star which then will be confirmed or not by follow-up observations[31].

Some examples of ground-based searches are the HATNet (Hungarian Automated Telescope Network) project which consists of seven telescopes located in three observatories and covers ∼ 10% of the sky, the OGLE (Optical Gravitational Lensing Experiment) the telescope of which was also used for photometric observations of stars in search of transiting planets, the SuperWASP (Wide Angle Search for Planets) survey which consists of two robotic observatories that cover both hemispheres, the KELT (Kilodegree Extremely Little Telescope) survey that uses two robotic telescopes searching for extrasolar planets around bright stars, while two of the most known space telescopes (past and current) that search for exoplanets are CoRoT (COnvection ROtation and planetary Transits) one of the two objectives of which was to look for transiting exoplanets and Kepler which was launched in 2009[31][17].

The importance of transit observations lies in the determination of the radii of the planets as well as of the masses when combined with radial velocity measurements. The stellar radius R _S , stellar mass M _S , planet radius R _p , semi-major axis a and inclination i are the five quantities which can be derived from a transit lightcurve since there is one unique solution. However, some assumptions must be made such as a circular planetary orbit, a known stellar mass-radius relation, M _p  M _S and the absence of an additional blended star[39].

A transit event is described by five equations. Three out of five equations illustrate the transit

depth, duration and shape. The two additional equations that are needed for the event to be completely

described are Kepler’s third law and the mass-radius relation of the star[39]. When a transit event is

observed, four quantities are documented, namely the period P, the transit depth ∆F , the duration

between the first and the fourth contacts t _T and the duration between the second and the third contacts

t _F , where the first, second, third and fourth contacts are shown in Fig. 9[39][31]. These four observables

are used in the three equations which describe the transit[39].

5.6 Transits 22

Figure 8: Depiction of a transit of a planet in front of a star when it is viewed edge-on[44].

∆F and t _T are given by[39]:

∆F = R _p R _S

! 2

(22)

t _T = P

π arcsin R _S a

( [1 + ( _R ^R

S

)] ² − [( _R ^a

S

) cos i] ² 1 − cos ² i

) 1/2 !

(23) The third equation that shows the transit shape is basically the ratio of the time interval of the flat part of the transit (t _F ) to the whole transit interval[39]:

sin (t _F π/P )

sin (t _T π/P ) = {[1 − (R _p /R _S )] ² − [(a/R _S ) cos i] ² } ^1/2

{[1 + (R _p /R _S )] ² − [(a/R _S ) cos i] ² } ^1/2 (24) The existence of 3458 as of 2 Mar 2017[17](3586 as of 7 Mar 2017[5]) extrasolar planets has been confirmed so far, 578[17] (603[5]) of which belong to multi-planet systems. The majority of these systems have been discovered either by Doppler measurements or by transits[14]. The study of planetary systems which consist of several planets reveals the dynamical stability and the interactions between the planets.

It also offers the chance for comparison to the Solar System and even for drawing conclusions about their properties from statistical data in the future.

In Fig.10 the transits of the six planets of the system Kepler-11 are shown. The inner five planets

compose a very densely packed system with their orbital periods ranging between 10 and 47 days. The

radii of the planets are 1.8 to 4.2 times the Earth radius[24].

6 EVOLUTION OF SOLAR-TYPE STARS 23

Figure 9: First, second, third and fourth contact of a transiting planet depicted by the numbers 1,2,3 and 4 respectively[39].

6 Evolution of solar-type stars

6.1 The H-R diagram

The initial composition of a star in terms of mass is 70% hydrogen, 25-30% helium and small quantities of heavier elements such as oxygen, carbon and nitrogen (CNO group)[11]. However, the evolution of each star depends on its initial mass. One tool that can help us study how different types of stars evolve is the H-R diagram, the two axes of which depict the surface temperature and the luminosity (or properties related to these quantities) of stars and can be deduced from observations[11].

In Fig.11 an H-R diagram is shown, where each point represents a star. Most of the stars follow

a line that runs diagonally from the left upper corner to the right bottom one of the diagram. These

stars are known as Main Sequence stars and they spend most of their lives in this status. Above the

Main Sequence and to the right, brighter stars are observed. If we draw a vertical line on the diagram,

for example at 4000K, all the stars on this line will have the same surface temperature. However, the

stars which are located above the Main Sequence will be brighter, despite the same value of surface

temperature. Furthermore, if we take into consideration the following equation which describes the

luminosity of a star with respect to its surface temperature and its radius:

6.1 The H-R diagram 24

Figure 10: Transits of the six planets b,c,d,e,f and g of the system Kepler-11[23].

6.2 The Main-Sequence (MS) 25

L = 4πR ² σT _{ef f} ² (25)

we see that the luminosity changes proportionally to radius. Thus, when two stars have the same effective temperature but one of them has a larger value of luminosity, the brighter one will also have a greater radius.

What is more, since they have low temperatures their spectra are shifted towards the red side of the spectrum, hence named red giants[11]. The lower left region of Fig.11, is populated by objects with high temperatures and low luminosities. These objects are called white dwarfs and constitute one of the final stages of stars. They are quite dense objects, since their radii is similar to the radius of Earth and their masses are comparable to the mass of the Sun[11].

Figure 11: Hertzsprung-Russell diagram[1]

6.2 The Main-Sequence (MS)

The range of different luminosities and temperatures we observe in the Main Sequence of the H-R diagram as well as the specific path that a star follows during its evolution depend on the initial mass of stars whereas the factor that determines if a star will be on the Main Sequence or not is its age[11].

Stars can be categorized in three groups based on their mass and thus their evolution, namely low-mass,

intermediate-mass and massive stars the mass of which ranges between the limits 0.7M  < M < 2M  ,

2M  < M < 9 − 10M  , M > 10M  respectively[11]. Stars with mass M < 0.7M  are known as red

dwarfs and even if their age is equal to the age of the universe they still are in the Main Sequence

6.2 The Main-Sequence (MS) 26

phase[11]. The lifetime of the Sun on the Main Sequence is 10 ¹⁰ yr while the Main Sequence lifetimes of a 0.5M  and of a 25M  stars are 7 × 10 ¹⁰ and 6 × 10 ⁶ years respectively[11].

The way a star counterbalances the gravitational pull of its own mass is by thermonuclear reactions inside its core. During the longest part of their lives, stars burn hydrogen in their core, which constitutes them as Main Sequence stars[11]. Through two chains of reactions, namely the proton-proton (p-p) chain and the CNO cycle, hydrogen vanishes gradually from the stellar core and is replaced mainly by helium[11]. The energy released by the fusion of four hydrogen nuclei to a helium nucleus, without taking into consideration the in-between steps, is equal to 26.731MeV[16]. p-p chain consists of three branches, p-p I, II and III[11][36].

p-p I

p + p → ² D + e ⁺ + ν (26)

2 D + p → ³ He + γ (27)

3 He + ³ He → ⁴ He + 2p (28)

p-p II

3 He + ⁴ He → ⁷ Be + γ (29)

7 Be + e ⁻ → ⁷ Li + ν (30)

7 Li + p → 2 ⁴ He (31)

The first reaction of p-p I, describes the encounter of two protons which results in deuterium through weak interactions. The second reaction shows the capture of a proton and the formation of ³ He. p-p chain I ends with the formation of a ⁴ He nucleus. In p-p II, the formation of ⁷ Be by two helium isotopes an electron capture, a proton capture and the production of two ⁴ He nuclei occur. Finally, p-p III, and consequently the p-p chain in total, is completed by the following reactions:

p-p III

3 He + ⁴ He → ⁷ Be + γ (32)

7 Be + p → ⁸ B + γ (33)

8 B → ⁸ Be + e ⁺ + ν (34)

8 Be → 2 ⁴ He (35)

6.3 The Red Giant Branch (RGB) 27

This last branch of the p-p chain describes the proton capture by ⁷ Be and the formation of ⁸ B, which decays to ⁸ Be and then breaks to two 2 ⁴ He nuclei. For temperatures between 5 × 10 ⁶ K < T < 8 × 10 ⁶ K the reaction rate of ³ He production is greater than the rate with which it is consumed, thus its abundance is high. When the temperature surpasses this upper limit though, the last reaction of p-p I and the first from p-p II and III become effective and its abundance decreases towards an equilibrium abundance[36].

Due to the presence even of small quantities of carbon, nitrogen and oxygen which act as catalysts, another process of hydrogen burning, the CNO cycle, occurs. There are two chains now, the CN cycle and the NO cycle[11][36]:

CN cycle

12 C + ¹ H → ¹³ N + γ

NO cycle

14 N + ¹ H → ¹⁵ O + γ

13 N → ¹³ C + e ⁺ + ν ¹⁵ O → ¹⁵ N + e ⁺ + ν

13 C + ¹ H → ¹⁴ N + γ ¹⁵ N + ¹ H → ¹⁶ O + γ

14 N + ¹ H → ¹⁵ O + γ ¹⁶ O + ¹ H → ¹⁷ F + γ

15 O → ¹⁵ N + e ⁺ + ν ¹⁷ F → ¹⁷ O + e ⁺ + ν

15 N + ¹ H → ¹² C + ⁴ He ¹⁷ O + ¹ H → ¹⁴ N + ⁴ He

The energy which is produced by thermonuclear reactions is transferred to outer layers either by radiation or by convection. The convective zone of stars with mass comparable to the Sun covers only a narrow shell below the photosphere approximately 2% of the star’s mass[11].

6.3 The Red Giant Branch (RGB)

When the hydrogen in a star’s core is depleted, the temperature is not high enough for helium burning to start. However, a thick hydrogen shell of ∼ 0.2M  , that surrounds the core has reached the neces- sary temperature for hydrogen burning[36]. The density inside the He core is high -in comparison to more massive stars at this state- and the core becomes electron degenerate. The degeneracy pressure counterbalances the gravitational pull of the overlying shells[36].

The transition from the MS to the RG phase passes through the Sub-Giant phase during which the stellar luminosity does not change significantly (Fig. 12) and the hydrogen shell becomes thinner due to H-exhaustion at the base of the shell[36]. When a star reaches the RGB, the shell width is of the order

∼ 0.001M  and it will decrease further down to ∼ 0.0001M  as the star reaches the tip of the RGB[36].

There is a strong correlation between the mass of the He core and the surface luminosity of a star as

it evolves through the RG phase which is described by the so-called He core mass-luminosity relation.

6.4 Mass-loss 28

The surface luminosity is determined by the properties of the hydrogen burning shell whose reaction rates depend on the the mass and radius of the degenerate core. Thus, an increase of the core’s mass will result in an greater value for the surface luminosity[36].

During the RG phase a star loses mass in a much greater rate than the previous stage, the Main Sequence, however, the core mass is not affected neither is the H-burning shell and consequently the luminosity of the star. The radius is affected by the mass-loss though and the size of the star is changing during RG phase[36]. As the star evolves its He core mass increases as does its temperature. The RG phase ends at the tip of the Red Giant Branch where the mass and the temperature of the core have reached the values M _He ∼ 0.48 − 0.50M  , T _c ∼ 10 ⁸ K respectively, causing the helium in the core to ignite. This results in a thermal runaway due to core instability caused by the degenerate electrons which is called He flash[36].

Figure 12: Evolutionary tracks of stars with different initial masses[3]

6.4 Mass-loss

Mass-loss occurs during all stages of stellar evolution. However, the rate differs from one stage to another

and it depends mainly on the stellar mass. Massive stars lose an important fraction of their total mass

on the Main Sequence, whereas mass-loss for low mass stars does not affect significantly their mass. The

7 SYMPLECTIC INTEGRATORS FOR N-BODY PROBLEMS 29

solar wind rate, for example, is of the order of ∼ 10 ⁻¹⁴ M  yr ⁻¹ [11].

Stars with mass ∼ 1M  spend ∼ 0.5Gyr in the RGB and when they reach the RGB tip, their mass has been reduced by ∼ 30%[32]. In order to investigate the effects of mass loss to the planetary orbits while a solar-type star is in the RGB we assume that the greater part of the mass is lost during the last million years in this stage. For the second part of this thesis, we implemented a constant mass-loss rate of 10 ⁻⁸ M  /yr for 30 million years to Kepler-11 and to the Sun.

7 Symplectic Integrators for N-Body Problems

The n-body problem, identified already by Newton, states that it is impossible to describe the motion of more than three objects, which interact gravitationally, analytically. Even for three bodies, analytical solutions exist only in special cases. Thus, numerical methods such as symplectic integrators are used to study the time evolution of planetary systems. The advantages of symplectic integrators over other ones are that the energy error is not built up over long term integrations and also that they are faster than other integrators when one of the bodies in the problem at hand is much more massive compared to the rest of the bodies, which is usually the case in planetary systems[9].

7.1 Hamiltonian systems

The dynamical evolution of planetary systems can be described by Hamilton’s equations. A Hamiltonian system is defined by the following differential equations:

dp i

dt = − ∂H

∂q _i (26)

dq _i

dt = ∂H

∂p i

(27) where H is the Hamiltonian function which must be real, sufficiently smooth and defined in a non- empty, open and connected set Y -the phase space- in the oriented Euclidean space R ^2d of (p,q) = (p ₁ , ...., p _d ; q ₁ , ...., q _d ), and p and q are the generalized momenta and coordinates respectively with the same dimension d which shows us the degrees of freedom of the system[37]. The function H usually describes the total mechanical energy of a system and when it can be written as a sum of kinetic and potential energy then it is called separable[37]. Hamiltonians can also be time-dependent (non- autonomous) H = H(p, q; t)[37].

The flow of a Hamiltonian system is defined as a smooth map: (H, t) → φ _t,H where t ∈ R. Let the initial conditions at t = 0 for a solution of the Hamiltonian equations be (p ⁰ , q ⁰ ) in Y. The value of this solution at time t is the flow of this set of initial conditions[37]:

(p,q) = φ _t,H (p ⁰ , q ⁰ ) (28)

When t is variable and (p ⁰ , q ⁰ ) is not, we return to the solution that was derived for the initial

conditions (p ⁰ , q ⁰ ). However, when we assume that (p ⁰ , q ⁰ ) varies but t is fixed, the flow (φ _t,H (p ⁰ , q ⁰ ))

can be defined only if the solution which corresponds to the initial conditions (p ⁰ , q ⁰ ), exists at t. The

7.2 Symplecticness 30

domain in which we define the flow will in general be smaller than Y, because the solutions might reach the boundary of Y in finite time, if |t| is large, and thus exist for bounded intervals of time[37].

7.2 Symplecticness

Numerical methods, such as integrators, are called symplectic when they are applied in Hamiltonian systems and the numerical solutions of them inherit the property of symplecticness[37]. In order to understand symplecticness, we start by studying it in systems with one degree of freedom in which it can be interpreted as preservation of area. The flow φ _t,H for real t, is an area-preserving transformation in Y[37]. This property of flow becomes clear through Liouville’s theorem which states that the phase space distribution function is constant along the trajectories of a system. The evolution of the phase space distribution ρ(p, q) is given by the Liouville equation:

dρ dt = ∂ρ

∂t + ∂ρ

∂q q + ˙ ∂ρ

∂p p ˙

!

= 0 (29)

where ˙ q and ˙ p are the time derivatives of q and p respectively. Assuming that the continuity equation can be applied to ρ(p, q), we get:

∂ρ

∂t + ∂(ρ ˙ q)

∂q + ∂(ρ ˙ p)

∂p

!

= 0 (30)

∂ρ

∂t + ∂ρ

∂q q + ˙ ∂ρ

∂p p ˙

!

+ ρ ∂ ˙ q

∂q + ∂ ˙ p

∂p

!

= 0

∂ρ

∂t + ∂ρ

∂q q + ˙ ∂ρ

∂p p ˙

!

+ ρ ∂

∂q

∂H

∂p − ∂

∂p

∂H

∂q

!

= 0 (31)

The terms ^∂ρ _∂t + 

∂ρ

∂q q + ˙ ^∂ρ _∂p p ˙ 

of Eq. 31 are the terms of Liouville equation. Thus, the Liouville equation is true when the last term of Eq. 31 is equal to 0:

∂

∂q

∂H

∂p − ∂

∂p

∂H

∂q

!

= 0 (32)

which means that the vector field [−∂H/∂q, ∂H/∂p] ^T is divergence free[37].

The next step is to check the area preservation, either by the Jacobian or by differential forms. Let ψ be a C ¹ transformation in a domain Y, (p ^∗ , q ^∗ ) = ψ(p, q). ψ is area preserving if and only if the Jacobian determinant is equal to 1, ∀(p, q) ∈ Y [37]:

∂p ^∗

∂p

∂q ^∗

∂q − ∂p ^∗

∂q

∂q ^∗

∂p = 1 (33)

If the Jacobian matrix of the transformation is ψ ⁰ = ^∂(p _∂(p,q)

^,q

⁾ , then we can write:

7.2 Symplecticness 31

ψ ^0T J ψ ⁰ = ∂(p ^∗ , q ^∗ )

∂(p, q)

T

J ∂(p ^∗ , q ^∗ )

∂(p, q) = J (34)

where J is the matrix:

J =  0 1

−1 0



Now, we construct a parallelogram P, using as a starting vertex the fixed point (p, q) and drawing two vectors v and w from this point. The parallelogram ψ(P) can be approximated by P ^∗ with ψ(p, q) being the fixed point and ψ ⁰ v, ψ ⁰ w the two vectors. The condition for two parallelograms, P and P ^∗ to have the same area is[37]:

v ^T ψ ^0T J ψ ⁰ w = v ^T J w (35)

If Eq.(34) is true then Eq.(35) is also true for all such parallelograms. Eq.(35) says that ∀(p, q) ∈ Y, the linear transformation ψ maps parallelograms with starting vertex the fixed point (p, q) to parallelo- grams with ψ(p, q) while the oriented area remains constant[37].

The other way to check the preservation of area in transformations is by using differential forms. In general, if P and Q are smooth real-valued functions and P (p, q)dp + Q(p, q)dq is a differential 1-form, the differential 1-forms dp ^∗ and dq ^∗ are[37]:

dp ^∗ = ∂p ^∗

∂p dp + ∂p ^∗

∂q dq (36)

dq ^∗ = ∂q ^∗

∂p dp + ∂q ^∗

∂q dq (37)

The exterior product of two differential 1-forms is a differential 2-form[37]:

dp ^∗ ∧ dq ^∗ = ∂p ^∗

∂p

∂q ^∗

∂p dp ∧ dp + ∂p ^∗

∂p

∂q ^∗

∂q dp ∧ dq + ∂p ^∗

∂q

∂q ^∗

∂p dq ∧ dp + ∂p ^∗

∂q

∂q ^∗

∂q dq ∧ dq (38) If we take into consideration that

dp ∧ dp = dq ∧ dq = 0 dp ∧ dq = −dq ∧ dp the differential 2-form becomes:

dp ^∗ ∧ dq ^∗ = ∂p ^∗

∂p

∂q ^∗

∂q − ∂p ^∗

∂q

∂q ^∗

∂p

!

dp ∧ dq (39)

If we recall Eq.(33), we see that a transformation ψ is area preserving if and only if the Jacobian determinant is equal to 1. Thus, if we apply this condition to Eq.(39), the relation that must be true for the preservation of area is[37]:

dp ^∗ ∧ dq ^∗ = dp ∧ dq (40)

The condition for area preservation in the form of Eq. (40) is usually more convenient in practice,

requiring fewer algebraic manipulations[37].