newtonian_energy_angular_momentum Subroutine

subroutine newtonian_energy_angular_momentum(eccentricity, periastron, mass1, mass2, energy, angular_momentum)

Compute the Newtonian energy and angular momentum of the system, imposing that the radial velocity of the fictitious body is 0 at the desired periastron, with the desired eccentricity.

The formulas used here are found by solving the equations that can be found in:

Goldstein, Poole, Safko, "Classical mechanics", Sec.3.2, eq.(3.16) with , and Sec.3.7, eq.(3.57) See Landau, Lifshitz, "Mechanics", Chapter III, eq.(14.5) with

for the energy and the angular momentum.

FT 16.12.2022

Arguments

Type IntentOptional Attributes Name
double precision, intent(in) :: eccentricity
double precision, intent(in) :: periastron
double precision, intent(in) :: mass1
double precision, intent(in) :: mass2
double precision, intent(out) :: energy
double precision, intent(out) :: angular_momentum

Called by

proc~~newtonian_energy_angular_momentum~~CalledByGraph proc~newtonian_energy_angular_momentum newtonian_energy_angular_momentum program~construct_newtonian_binary construct_newtonian_binary program~construct_newtonian_binary->proc~newtonian_energy_angular_momentum

Contents

Variables

mu

Source Code

newtonian_energy_angular_momentum

Variables

Type Visibility Attributes Name Initial
double precision, public :: mu

Source Code

  SUBROUTINE newtonian_energy_angular_momentum &
    (eccentricity, periastron, mass1, mass2, energy, angular_momentum)

    !***********************************************************
    !
    !# Compute the Newtonian energy and angular momentum of the system,
    !  imposing that the radial velocity of the fictitious
    !  body is 0 at the desired periastron, with the desired eccentricity.
    !
    !  The formulas used here are found by solving the equations
    !  that can be found in:
    !
    !  Goldstein, Poole, Safko, "Classical mechanics",
    !  Sec.3.2, eq.(3.16) with \(\dot(r)=0\), and Sec.3.7, eq.(3.57)
    !  See Landau, Lifshitz, "Mechanics", Chapter III,
    !  eq.(14.5) with \(\dot(r)=0\)
    !
    !  for the energy and the angular momentum.
    !
    !  FT 16.12.2022
    !
    !***********************************************************

    IMPLICIT NONE

    DOUBLE PRECISION, INTENT(IN) :: eccentricity, periastron, mass1, mass2

    DOUBLE PRECISION, INTENT(OUT):: energy, angular_momentum

    DOUBLE PRECISION:: mu

    mu= mass1*mass2/(mass1 + mass2)

    IF(eccentricity == zero)THEN
    ! Circle

      angular_momentum= SQRT(mu*mass1*mass2*periastron)

    ELSEIF(eccentricity == one)THEN
    ! Parabola (straight line is not considered here; it would have zero
    ! angular momentum)

      angular_momentum= SQRT(two*mu*mass1*mass2*periastron)

    ELSEIF(eccentricity > one)THEN
    ! Hyperbola

      angular_momentum= SQRT((one + eccentricity)*mu*mass1*mass2*periastron)

    ELSEIF(zero < eccentricity .AND. eccentricity < one)THEN
    ! Ellipse [SQRT((one - eccentricity)*mu*mass1*mass2*periastron) would be
    ! for an ellispe having apoastron equal to our value of the periastron]

      angular_momentum= SQRT((one + eccentricity)*mu*mass1*mass2*periastron)

    ELSE

      PRINT *, "** ERROR in SUBROUTINE newtonian_energy_angular_momentum!"
      PRINT *, " * The value for the eccentricity is negative!"
      PRINT *, "   eccentricity= ", eccentricity
      PRINT *, " * Stopping..."
      PRINT *
      STOP

    ENDIF

    energy= &
      mu*(mass1*mass2)**2*(eccentricity**2 - one)/(two*angular_momentum**2)

  END SUBROUTINE newtonian_energy_angular_momentum