Commit 3c50535e authored by Romolo Politi's avatar Romolo Politi
Browse files

Upload Version 1.0

parents
# DPI Library - Dynamical Plug-In for Mercury 6 (release 1.0)
## Author: Diego Turrini
DPI is a FORTRAN77 library that supplies the symplectic mapping method for binary star systems (Chambers, Quintana, Duncan & Lissauer, 2002, AJ 123, 2884-2894) for the Mercury N-Body software package by John E. Chambers (ascl:1201.008). The binary symplectic mapping is implemented as a hybrid symplectic method (Chambers, 1999, MNRAS 304, 793-799) that allows close encounters and collisions between massive bodies and is therefore suitable for planetary accretion simulations.
License: DPI is released under GPL 3 license.
Note on this release: release 1.0 supplies the symplectic mapping for S-type binary systems (also called wide binary systems in Chambers et al. 2002), where the planets orbit one of the stars and are perturbed by the stellar companion (which is the outermost body of the system).
Disclaimer: the DPI library is supplied with a modified version of Mercury 6 for the convenience of the users; however, only those files whose name begins with "dpi" are property of the author and are released under GPL 3 license.
## USING DPI LIBRARY WITH MERCURY 6
The use of DPI in conjunction with Mercury 6 should be straightforward for Mercury's users. Users not familiar with Mercury 6 are advised to first read Mercury's manual to get in touch with its working.
In order to use Mercury together with DPI the users need to compile together the following two files:
- dpi1_0.for, containing the DPI library
- mercury6_2-dpi.for, containing a modified version of Mercury 6 adapted to call the subroutines provided by the DPI library
The two files can either be compiled together, e.g:
gfortran mercury6_2-dpi.for dpi-rc2.for -o executablename
or separately by using the option "-c" to create binary object files, e.g.:
gfortran -c mercury6_2-dpi.for
gfortran -c dpi1_0.for
gfortran mercury6_2-dpi.o dpi1_0.o -o executablename
Once compiled, to use the S-type symplectic mapping:
- create a file big.in with the bodies you want to integrate, taking care that the binary companion is the last body listed
- indicate "wb" (wide binary, since the binary star is the outermost body of the system) as the algorithm in the file param.in of Mercury
If the program is run using the test input files supplied with DPI, the energy error should oscillate around 1E-7 and the angular momentum error varying between 1E-13 and 1E-14 (the behaviour of the latter can vary depending on the hardware and the compiler due to the actual precision of the computations and the way rounding errors are handled).
)O+_05 Big-body initial data
) Lines beginning with ) are ignored
)--------------------------------------------------
style (Cartesian,Asteroidal,Cometary) = Asteroidal
epoch (in days) = 0.0
)--------------------------------------------------
Planet1 m= 1.0d-3 r= 1. d= 1.33000004
5.000E+00 0.000E-00 0.000E+00
0.000E+00 0.000E+00 0.000E+00
0.000E+00 0.000E+00 0.000E+00
Planet2 m= 1.0d-3 r= 1. d= 1.33000004
1.000E+01 0.000E-00 0.000E+00
0.000E+00 0.000E+00 0.000E+00
0.000E+00 0.000E+00 0.000E+00
Planet3 m= 1.0d-3 r= 1. d= 1.33000004
2.000E+01 0.000E-00 0.000E+00
0.000E+00 0.000E+00 0.000E+00
0.000E+00 0.000E+00 0.000E+00
\ No newline at end of file
This diff is collapsed.
)O+_06 element (WARNING: Do not delete this line!!)
) Lines beginning with `)' are ignored.
)---------------------------------------------------------------------
number of input files = 1
)---------------------------------------------------------------------
) List the input files, one per line
xv.out
)---------------------------------------------------------------------
type of elements (central body, barycentric, Jacobi) = ce
minimum interval between outputs (days) = 100.d0
express time in days or years = years
express time relative to integration start time = yes
)---------------------------------------------------------------------
) Output format? (e.g. a8.4 => semi-major axis with 8 digits & 4 dec. places)
a21.16 e21.16 i21.16 x21.16 y21.16 z21.16 u21.16 v21.16 w21.16
)---------------------------------------------------------------------
) Which bodies do you want? (List one per line or leave blank for all bodies)
) a21.16 e8.6 i8.4 g8.4 n8.4 l8.4 m13e
This diff is collapsed.
big.in
small.in
param.in
xv.out
ce.out
info.out
big.dmp
small.dmp
param.dmp
restart.dmp
c%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
c
c MERCURY.INC (ErikSoft 4 March 2001)
c
c%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
c
c Author: John E. Chambers
c
c Parameters that you may want to alter at some point:
c
c NMAX = maximum number of bodies
c CMAX = maximum number of close-encounter minima monitored simultaneously
c NMESS = maximum number of messages in message.in
c HUGE = an implausibly large number
c NFILES = maximum number of files that can be open at the same time
c
integer NMAX, CMAX, NMESS, NFILES
real*8 HUGE
c
parameter (NMAX = 1000)
parameter (CMAX = 1000)
parameter (NMESS = 200)
parameter (HUGE = 9.9d29)
parameter (NFILES = 90)
c
c------------------------------------------------------------------------------
c
c Constants:
c
c DR = conversion factor from degrees to radians
c K2 = Gaussian gravitational constant squared
c AU = astronomical unit in cm
c MSUN = mass of the Sun in g
c
real*8 PI,TWOPI,PIBY2,DR,K2,AU,MSUN
c
parameter (PI = 3.141592653589793d0)
parameter (TWOPI = PI * 2.d0)
parameter (PIBY2 = PI * .5d0)
parameter (DR = PI / 180.d0)
parameter (K2 = 2.959122082855911d-4)
parameter (AU = 1.4959787d13)
parameter (MSUN = 1.9891d33)
This diff is collapsed.
1 6 days
2 6 years
3 13 solar masses
4 3 AU
5 3 no
6 3 yes
7 3 low
8 6 medium
9 4 high
10 0
11 33 Integration parameters
12 33 ----------------------
13 14 Algorithm:
14 38 Second-order mixed-variable symplectic
15 24 Bulirsch-Stoer (general)
16 37 Bulirsch-Stoer (conservative systems)
17 16 15th-order RADAU
18 0
19 0
20 0
21 0
22 5 Test
23 48 Hybrid symplectic integrator (mixed coordinates)
24 44 Hybrid symplectic (close binary coordinates)
25 43 Hybrid symplectic (wide binary coordinates)
26 32 Integration start epoch:
27 32 Integration stop epoch:
28 32 Output interval:
29 32 Element origin:
30 31 Initial timestep:
31 36 Accuracy parameter:
32 36 Central mass:
33 36 J_2:
34 36 J_4:
35 36 J_6:
36 36 Ejection distance:
37 36 Radius of central body:
38 29 Number of Big bodies:
39 29 Number of Small bodies:
40 37 Output precision:
41 40 Includes collisions:
42 40 Includes fragmentation:
43 0
44 0
45 40 Includes relativity:
46 40 Includes user-defined force routine:
47 10 barycentre
48 12 central body
49 0
50 0
51 30 Integration details
52 30 -------------------
53 29 Initial energy:
54 29 Initial angular momentum:
55 65 Integrating massive bodies and particles up to the same epoch.
56 34 Beginning the main integration.
57 24 Integration complete.
58 48 Fractional energy change due to integrator:
59 48 Fractional angular momentum change:
60 57 Fractional energy change due to collisions/ejections:
61 57 Fractional angular momentum change:
62 47 Continuing integration from dump files at
63 6 Time:
64 6 Date:
65 9 dE/E:
66 9 dL/L:
67 35 collided with the central body at
68 12 ejected at
69 12 was hit by
70 34 removed due to an encounter with
71 4 at
72 26 solar masses AU^2 day^-2
73 26 solar masses AU^2 day^-1
74 36 lost mass due to rotational breakup
75 24 removed due to small a
76 0
77 0
78 0
79 0
80 0
81 8 ERROR:
82 49 Modify mercury.inc and recompile Mercury.
83 62 Check the file containing initial data for Big bodies.
84 64 Check the file containing initial data for Small bodies.
85 57 Check the file containing integration parameters.
86 22 Check files.in
87 27 This file already exists:
88 34 This file is needed to continue:
89 30 This filename is duplicated:
90 40 The total number of bodies exceeds NMAX.
91 68 Data style on first line must be Cartesian, Asteroidal or Cometary
92 68 You cannot integrate non-gravitational forces using this algorithm.
93 64 You cannot integrate a user-defined force using this algorithm.
94 64 You cannot integrate massive Small bodies using this algorithm.
95 66 Massive Small bodies must have the same epoch as the Big bodies.
96 49 Check character implies input file is corrupted.
97 62 Mass, density, encounter limit must be >= 0 for this object:
98 46 This integration algorithm is not available:
99 50 A problem occurred reading the parameter on line
100 50 A problem occurred reading data for this object:
101 56 A problem occured reading the epoch for the Big bodies.
102 67 You cannot use non-zero J2,J4,J6 using the close-binary algorithm.
103 34 Two objects both have this name:
104 36 is corrupted at line number:
105 42 Central-body radius exceeds maximum radius.
106 68 Maximum/Central radius is large. Output precision will be degraded.
107 58 Coordinate origin must be Central, Barycentric or Jacobi.
108 0
109 0
110 0
111 0
112 0
113 0
114 0
115 0
116 0
117 0
118 0
119 0
120 0
121 10 WARNING:
122 53 Truncating the name of this object to 8 characters:
123 30 Main integration is backwards.
124 26 No Big bodies are present.
125 28 No Small bodies are present.
126 50 Stopping integration due to an encounter between
127 45 Throwing this object into the central body:
128 42 Setting output threshhold DA to infinity.
129 42 Setting output threshhold DE to infinity.
130 42 Setting output threshhold DI to infinity.
131 43 Increasing the radius of the central body.
132 56 Total number of current close encounters exceeds CMAX.
133 0
134 0
135 0
136 0
137 0
138 0
139 0
140 0
141 0
142 0
143 0
144 0
145 0
146 0
147 0
148 0
149 0
150 0
151 67 )O+_05 Integration parameters (WARNING: Do not delete this line!!)
152 66 )O+_05 Big-body initial data (WARNING: Do not delete this line!!)
153 68 )O+_05 Small-body initial data (WARNING: Do not delete this line!!)
154 39 ) Lines beginning with `)' are ignored.
155 70 )---------------------------------------------------------------------
156 43 style (Cartesian, Asteroidal, Cometary) =
157 20 epoch (in days) =
158 35 ) Important integration parameters:
159 48 algorithm (MVS, BS, BS2, RADAU, HYBRID etc) =
160 21 start time (days) =
161 20 stop time (days) =
162 26 output interval (days) =
163 19 timestep (days) =
164 22 accuracy parameter =
165 22 ) Integration options:
166 44 stop integration after a close encounter =
167 29 allow collisions to occur =
168 37 include collisional fragmentation =
169 33 express time in days or years =
170 51 express time relative to integration start time =
171 20 output precision =
172 24 < Not used at present >
173 37 include relativity in integration =
174 30 include user-defined force =
175 52 ) These parameters do not need to be adjusted often:
176 26 ejection distance (AU) =
177 31 radius of central body (AU) =
178 31 central mass (solar masses) =
179 14 central J2 =
180 14 central J4 =
181 14 central J6 =
182 24 < Not used at present >
183 24 < Not used at present >
184 45 Hybrid integrator changeover (Hill radii) =
185 42 number of timesteps between data dumps =
186 48 number of timesteps between periodic effects =
187 41 origin (Central, Barycentric, Jacobi) =
188 0
189 0
190 0
191 0
192 0
193 0
194 0
195 0
196 0
197 0
198 0
199 0
200 0
)O+_05 Integration parameters (WARNING: Do not delete this line!!)
) Lines beginning with `)' are ignored.
)---------------------------------------------------------------------
) Important integration parameters:
)---------------------------------------------------------------------
algorithm (MVS, BS, BS2, RADAU, HYBRID etc) = wb
start time (days) = 0.
stop time (days) = 36525000
output interval (days) = 36525
timestep (days) = 50
accuracy parameter = 1.E-14
)---------------------------------------------------------------------
) Integration options:
)---------------------------------------------------------------------
stop integration after a close encounter = no
allow collisions to occur = yes
include collisional fragmentation = no
express time in days or years = years
express time relative to integration start time = yes
output precision = high
output barycentric instead of normal elements = no
include relativity in integration = no
include user-defined force = no
)---------------------------------------------------------------------
) These parameters do not need to be adjusted often:
)---------------------------------------------------------------------
ejection distance (AU) = 5000.
radius of central body (AU) = 0.0007
central mass (solar masses) = 1.0
central J2 = 0.
central J4 = 0.
central J6 = 0.
maximum eccentricity used in output = 1.1
minimum pericentric distance (units of rcen) = 0.1
Hybrid integrator changeover (Hill radii) = 3.0
number of timesteps between data dumps = 36525
number of timesteps between periodic effects = 36525
)O+_06 Small-body initial data (WARNING: Do not delete this line!!)
) Lines beginning with `)' are ignored.
)---------------------------------------------------------------------
style (Cartesian, Asteroidal, Cometary) = Ast
)---------------------------------------------------------------------
Test
1.000E+00 0.000E-00 0.000E+00
0.000E+00 0.000E+00 0.000E+00
0.000E+00 0.000E+00 0.000E+00
\ No newline at end of file
c*************************************************************************
c SWIFT.INC
c*************************************************************************
C Include file for SWIFT
c
c Author: Hal Levison
c Date: 2/2/93
c Last revision: 3/7/93
implicit NONE ! you got it baby
c... Maximum array size
integer NPLMAX, NTPMAX
parameter (NPLMAX = 202) ! max number of planets, including the Sun
parameter (NTPMAX = 2000) ! max number of test particles
c... Size of the test particle status flag
integer NSTAT
parameter (NSTAT = 3)
c... convergence criteria for danby
real*8 DANBYAC , DANBYB
parameter (DANBYAC= 1.0d-14, DANBYB = 1.0d-13)
c... loop limits in the Laguerre attempts
integer NLAG1, NLAG2
parameter(NLAG1 = 50, NLAG2 = 400)
c... A small number
real*8 TINY
PARAMETER(TINY=4.D-15)
c... trig stuff
real*8 PI,TWOPI,PIBY2,DEGRAD
parameter (PI = 3.141592653589793D0)
parameter (TWOPI = 2.0D0 * PI)
parameter (PIBY2 = PI/2.0D0)
parameter (DEGRAD = 180.0D0 / PI)
c-------------------------------------------------------------------------
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment