#! /usr/bin/env python

# planetaryOrbits.py is program which calculates the properties of an 
# eliptical planetary orbit.
#
# The program asks the user to input the distance to the Sun and the velocity 
# at the perhelion and then claculates and prints out the planetary distance 
# and velocity at the aphelion (most distant and slowest), the orbital period, 
# and eccentricity.  
#
#
# conserving energy
# 1/2 m v1**2 - GmM/l1 == 1/2 m v2**2 -GmM/l2
# v1**2 - 2GM/l1 == v2**2 - 2GM/l2
#    with l2 v2 == l1 v1
# v2**2 -(2GM/(l1*v1))*v2 -(v1**2 - 2GM/l1) == 0
#
# Using the quadratic formula we need the smaller solution
# as that is the one for aphelion(slowest velocity) the other
# solution would be for the perihelion velocity
# A == 1
# B == -2*G*M/(v1*l1)
# C == -v1**2 + 2*G*M/l1
#
#

### OUTPUT ###
#hpc-login-24 555% planetaryOrbits.py
#Enter the distance at the perihelion in meters: 1.471e11
#Enter the velocity at the perihelion in meters per second: 3.0287e4
#
#Distance at the aphelion is 152027197209 m
#Velocity at the aphelion is 29305 m/s.
#Period: 31543060 seconds, 365 days,  1 years.
#Eccentricity: 0.0165
#
#hpc-login-24 554% planetaryOrbits.py
#Enter the distance at the perihelion in meters: 8.783e10
#Enter the velocity at the perihelion in meters per second: 5.4529e4
#
#Distance at the aphelion is 5282214660876 m
#Velocity at the aphelion is 907 m/s.
#Period: 2399312512 seconds, 27770 days, 76 years.
#Eccentricity: 0.9673


# Paul Eugenio 
# PHZ4151C
# Jan 18, 2018

# program header code
from __future__ import division, print_function
from math import pi, sqrt

# main body of program
G = 6.6738e-11     # Newton's gravitional constant
M = 1.9891e30      # Mass of the Sun (kg)

l1 = float(raw_input("Enter the distance at the perihelion in meters: "))
v1 = float(raw_input("Enter the velocity at the perihelion in meters per second: "))

# calculate the distance at the aphelion
# A = 1 
B = -2*G*M/(v1*l1)
C = -v1**2 + 2*G*M/l1
v2 = (-B - sqrt(B*B - 4*C))/2
l2 = l1*v1/v2

# Calculate the other orbital parameters

a = (l1+l2)/2          # Semi-major axis
b = sqrt(l1*l2)        # Semi-minor axis
T = 2*pi*a*b/(l1*v1)   # Orbital period
e = (l2-l1)/(l2+l1)    # Orbital eccentricity

# Print out results
print()
print( "Distance at the aphelion is {0:2.0f} m".format(l2) )
print( "Velocity at the aphelion is {0:2.0f} m/s.".format(v2) )
print( "Period: {0:2.0f} seconds, {1:2.0f} days, {2:2.0f} years.".format(T, T/(24*3600), T/(24*3600*365.25)) ) 

print( "Eccentricity: {0:2.4f}".format(e) )