Nbody6
 All Files Functions Variables
circ.f
Go to the documentation of this file.
1  SUBROUTINE circ
2 *
3 *
4 * Circularization of inner binary.
5 * --------------------------------
6 *
7  include 'common6.h'
8  common/koztmp/ emax
9 *
10 *
11 * Determine inner eccentricity after merger termination.
12  ip = npairs
13  i1 = 2*ip - 1
14  a0 = -0.5*body(n+ip)/h(ip)
15  ecc2 = (1.0 - r(ip)/a0)**2 + tdot2(ip)**2/(body(n+ip)*a0)
16  ecc0 = sqrt(ecc2)
17 *
18 * Adopt circularization based on constant angular momentum.
19 * SEMI = A0*(1.0 - ECC2)
20  semi = a0*(1.0 - emax**2)
21  semi0 = semi
22  rx = max(radius(i1),radius(2*ip))
23 * Impose a limit to avoid Roche case (but note BRAKE2 reduction).
24  semi = max(semi,3.0*rx)
25  zmu = body(i1)*body(2*ip)/body(n+ip)
26  hi = h(ip)
27 * Set new binding energy and update ECOLL for conservation.
28  h(ip) = -0.5*body(n+ip)/semi
29  ecoll = ecoll + zmu*(hi - h(ip))
30 *
31 * Transform to new semi-major axis and small eccentricity (KSTAR=10).
32  CALL expand(ip,a0)
33  CALL ksperi(ip)
34  ecc1 = 0.001
35  CALL deform(ip,ecc0,ecc1)
36  kstar(n+ip) = max(10,kstar(n+ip))
37  tev0(n+ip) = time
38  tev(n+ip) = time + 0.1 ! important for exceeding TEV0.
39 *
40 * Estimate circularization time in Myr (C. Tout, Cambody lecture 2006).
41  IF (radius(2*ip).GT.radius(i1)) THEN
42  j1 = i1 + 1
43  j2 = i1
44  ELSE
45  j1 = i1
46  j2 = i1 + 1
47  END IF
48  q1 = body(j1)/body(j2)
49  tcirc = 2.0*q1**2/(1.0 + q1)*(semi/radius(j1))**8
50  tcirc = 1.0d-06*tcirc
51 *
52  WRITE (6,15) time+toff, name(i1), ecc0, a0, semi, tcirc
53  WRITE (45,15) time+toff, name(i1), ecc0, a0, semi, tcirc
54  15 FORMAT (' CIRCULARIZED T NM E A0 A TC ',f9.2,i6,f8.4,1p,3e10.2)
55  CALL flush(45)
56 *
57 * Check standard collision condition.
58  IF (radius(i1) + radius(2*ip).GT.semi) THEN
59  kspair = npairs
60  WRITE (6,20) kspair, ecc0, a0, semi, rx
61  20 FORMAT (' KOZAI COAL/COLL KS E0 A0 A RX ',
62  & i4,f8.4,1p,3e10.2)
63  iqcoll = -2
64  CALL cmbody(r(ip),2)
65  END IF
66 *
67  RETURN
68 *
69  END