anal_jacobi.f 2.0 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566
  1. c*************************************************************************
  2. c ANAL_JACOBI.F
  3. c*************************************************************************
  4. c Calculates the Jacobi integral of a particle given the position and velocity
  5. c and ang. velocity vector of the Sun and Planet
  6. c
  7. c Input:
  8. c msun ==> Mass of the Sun (real scalar)
  9. c mpl ==> Mass of the Planet (real scalar)
  10. c omegax,..y,z ==> vector components of planet-sun orbit's
  11. c angular velocity vector(real scalars)
  12. c xb,yb,zb ==> Barycentric position of particle
  13. c (real scalars)
  14. c vxb,vyb,vzb ==> Barycentric vel of particle (real scalars)
  15. c xsb,ysb,zsb ==> Barycentric position of Sun (real scalars)
  16. c xplb,yplb,zplb ==> Barycentric position of planet (real scalars)
  17. c
  18. c Output:
  19. c jacobi ==> Value of the jacobi constant
  20. c
  21. c Remarks:
  22. c Authors: Martin Duncan
  23. c Date: April 13/93
  24. c Last revision: 3/4/93 HFL
  25. subroutine anal_jacobi(msun,mpl,omegax,omegay,omegaz,xb,yb,zb,
  26. & vxb,vyb,vzb,xsb,ysb,zsb,xplb,yplb,zplb,jacobi)
  27. include '../swift.inc'
  28. c... Inputs:
  29. real*8 msun,mpl,omegax,omegay,omegaz,xb,yb,zb,vxb,vyb,vzb
  30. real*8 xsb,ysb,zsb,xplb,yplb,zplb
  31. c... Outputs:
  32. real*8 jacobi
  33. c... Internals
  34. real*8 rr,jx,jy,jz
  35. c----
  36. c... Executable code
  37. jacobi = 0.5d0*(vxb**2 + vyb**2 + vzb**2)
  38. rr = sqrt((xb-xsb)**2 + (yb-ysb)**2 + (zb-zsb)**2)
  39. jacobi = jacobi - msun/rr
  40. rr = sqrt((xb-xplb)**2 + (yb-yplb)**2 + (zb-zplb)**2)
  41. jacobi = jacobi - mpl/rr
  42. jx = yb*vzb - zb*vyb
  43. jy = zb*vxb - xb*vzb
  44. jz = xb*vyb - yb*vxb
  45. jacobi = jacobi - omegax*jx - omegay*jy - omegaz*jz
  46. return
  47. end ! anal_jacobi
  48. c-------------------------------------------------------------------------