src/amuse/ext/orbital_elements.py
def rel_posvel_arrays_from_orbital_elements(
        primary_mass,
        secondary_mass,
        semi_major_axis,
        eccentricity=0 | units.rad,
        true_anomaly=0 | units.rad,
        inclination=0 | units.rad,
        longitude_of_the_ascending_node=0 | units.rad,
        argument_of_periapsis=0 | units.rad,
        G=None
        ):
    """
    Returns relative positions/velocities for secondaries orbiting primaries.
    If primary_mass is a scalar, assumes the same primary for all secondaries.
    """
    if G is None:
        G=derive_G(primary_mass)

    try:
        number_of_secondaries = len(secondary_mass)
    except:
        number_of_secondaries = 1

    # arrays need to be equal to number of secondaries, or have just one value
    primary_mass = equal_length_array_or_scalar(
            primary_mass, length=number_of_secondaries)
    semi_major_axis = equal_length_array_or_scalar(
            semi_major_axis, length=number_of_secondaries)
    eccentricity = equal_length_array_or_scalar(
            eccentricity, length=number_of_secondaries)
    true_anomaly = equal_length_array_or_scalar(
            true_anomaly, length=number_of_secondaries)
    inclination = equal_length_array_or_scalar(
            inclination, length=number_of_secondaries)
    longitude_of_the_ascending_node = equal_length_array_or_scalar(
            longitude_of_the_ascending_node, length=number_of_secondaries)
    argument_of_periapsis = equal_length_array_or_scalar(
            argument_of_periapsis, length=number_of_secondaries)

    cos_true_anomaly = cos(true_anomaly)
    sin_true_anomaly = sin(true_anomaly)

    cos_inclination = cos(inclination)
    sin_inclination = sin(inclination)

    cos_arg_per = cos(argument_of_periapsis)
    sin_arg_per = sin(argument_of_periapsis)

    cos_long_asc_nodes = cos(longitude_of_the_ascending_node)
    sin_long_asc_nodes = sin(longitude_of_the_ascending_node)

    # alpha is a unit vector directed along the line of node
    alphax = (
            cos_long_asc_nodes*cos_arg_per
            - sin_long_asc_nodes*sin_arg_per*cos_inclination
            )
    alphay = (
            sin_long_asc_nodes*cos_arg_per
            + cos_long_asc_nodes*sin_arg_per*cos_inclination
            )
    alphaz = sin_arg_per*sin_inclination
    alpha = numpy.array([alphax, alphay, alphaz])

    # beta is a unit vector perpendicular to alpha and the orbital angular
    # momentum vector
    betax = (
            - cos_long_asc_nodes*sin_arg_per
            - sin_long_asc_nodes*cos_arg_per*cos_inclination
            )
    betay = (
            - sin_long_asc_nodes*sin_arg_per
            + cos_long_asc_nodes*cos_arg_per*cos_inclination
            )
    betaz = cos_arg_per*sin_inclination
    beta = numpy.array([betax, betay, betaz])

    # Relative position and velocity
    separation = (  # Compute the relative separation
            semi_major_axis*(1.0 - eccentricity**2)
            / (1.0 + eccentricity*cos_true_anomaly)
            )
    position_vector = (
            separation*cos_true_anomaly*alpha
            + separation*sin_true_anomaly*beta
            ).T
    velocity_tilde = (
            (
                G*(primary_mass + secondary_mass)
                / (semi_major_axis*(1.0 - eccentricity**2))
                )**0.5
            )  # Common factor
    velocity_vector = (
            -1.0 * velocity_tilde * sin_true_anomaly * alpha
            + velocity_tilde*(eccentricity + cos_true_anomaly)*beta
            ).T

    return position_vector, velocity_vector

    Returns relative positions/velocities for secondaries orbiting primaries.
    If primary_mass is a scalar, assumes the same primary for all secondaries.
    
Returns relative positions/velocities for secondaries orbiting primaries.
If primary_mass is a scalar, assumes the same primary for all secondaries.
Returns relative positions / velocities for secondaries orbiting primaries . If primary_mass is a scalar assumes the same primary for all secondaries .
def rel_posvel_arrays_from_orbital_elements(
        primary_mass,
        secondary_mass,
        semi_major_axis,
        eccentricity=0 | units.rad,
        true_anomaly=0 | units.rad,
        inclination=0 | units.rad,
        longitude_of_the_ascending_node=0 | units.rad,
        argument_of_periapsis=0 | units.rad,
        G=None
        ):
    if G is None:
        G=derive_G(primary_mass)
    try:
        number_of_secondaries = len(secondary_mass)
    except:
        number_of_secondaries = 1
    primary_mass = equal_length_array_or_scalar(
            primary_mass, length=number_of_secondaries)
    semi_major_axis = equal_length_array_or_scalar(
            semi_major_axis, length=number_of_secondaries)
    eccentricity = equal_length_array_or_scalar(
            eccentricity, length=number_of_secondaries)
    true_anomaly = equal_length_array_or_scalar(
            true_anomaly, length=number_of_secondaries)
    inclination = equal_length_array_or_scalar(
            inclination, length=number_of_secondaries)
    longitude_of_the_ascending_node = equal_length_array_or_scalar(
            longitude_of_the_ascending_node, length=number_of_secondaries)
    argument_of_periapsis = equal_length_array_or_scalar(
            argument_of_periapsis, length=number_of_secondaries)
    cos_true_anomaly = cos(true_anomaly)
    sin_true_anomaly = sin(true_anomaly)
    cos_inclination = cos(inclination)
    sin_inclination = sin(inclination)
    cos_arg_per = cos(argument_of_periapsis)
    sin_arg_per = sin(argument_of_periapsis)
    cos_long_asc_nodes = cos(longitude_of_the_ascending_node)
    sin_long_asc_nodes = sin(longitude_of_the_ascending_node)
    alphax = (
            cos_long_asc_nodes*cos_arg_per
            - sin_long_asc_nodes*sin_arg_per*cos_inclination
            )
    alphay = (
            sin_long_asc_nodes*cos_arg_per
            + cos_long_asc_nodes*sin_arg_per*cos_inclination
            )
    alphaz = sin_arg_per*sin_inclination
    alpha = numpy.array([alphax, alphay, alphaz])
    betax = (
            - cos_long_asc_nodes*sin_arg_per
            - sin_long_asc_nodes*cos_arg_per*cos_inclination
            )
    betay = (
            - sin_long_asc_nodes*sin_arg_per
            + cos_long_asc_nodes*cos_arg_per*cos_inclination
            )
    betaz = cos_arg_per*sin_inclination
    beta = numpy.array([betax, betay, betaz])
    separation = (  
            semi_major_axis*(1.0 - eccentricity**2)
            / (1.0 + eccentricity*cos_true_anomaly)
            )
    position_vector = (
            separation*cos_true_anomaly*alpha
            + separation*sin_true_anomaly*beta
            ).T
    velocity_tilde = (
            (
                G*(primary_mass + secondary_mass)
                / (semi_major_axis*(1.0 - eccentricity**2))
                )**0.5
            )  
    velocity_vector = (
            -1.0 * velocity_tilde * sin_true_anomaly * alpha
            + velocity_tilde*(eccentricity + cos_true_anomaly)*beta
            ).T
    return position_vector, velocity_vector
def rel_posvel_arrays_from_orbital_elements ( primary_mass , secondary_mass , semi_major_axis , eccentricity = 0 | units . rad , true_anomaly = 0 | units . rad , inclination = 0 | units . rad , longitude_of_the_ascending_node = 0 | units . rad , argument_of_periapsis = 0 | units . rad , G = None ) : if G is None : G = derive_G ( primary_mass ) try : number_of_secondaries = len ( secondary_mass ) except : number_of_secondaries = 1 primary_mass = equal_length_array_or_scalar ( primary_mass , length = number_of_secondaries ) semi_major_axis = equal_length_array_or_scalar ( semi_major_axis , length = number_of_secondaries ) eccentricity = equal_length_array_or_scalar ( eccentricity , length = number_of_secondaries ) true_anomaly = equal_length_array_or_scalar ( true_anomaly , length = number_of_secondaries ) inclination = equal_length_array_or_scalar ( inclination , length = number_of_secondaries ) longitude_of_the_ascending_node = equal_length_array_or_scalar ( longitude_of_the_ascending_node , length = number_of_secondaries ) argument_of_periapsis = equal_length_array_or_scalar ( argument_of_periapsis , length = number_of_secondaries ) cos_true_anomaly = cos ( true_anomaly ) sin_true_anomaly = sin ( true_anomaly ) cos_inclination = cos ( inclination ) sin_inclination = sin ( inclination ) cos_arg_per = cos ( argument_of_periapsis ) sin_arg_per = sin ( argument_of_periapsis ) cos_long_asc_nodes = cos ( longitude_of_the_ascending_node ) sin_long_asc_nodes = sin ( longitude_of_the_ascending_node ) alphax = ( cos_long_asc_nodes * cos_arg_per - sin_long_asc_nodes * sin_arg_per * cos_inclination ) alphay = ( sin_long_asc_nodes * cos_arg_per + cos_long_asc_nodes * sin_arg_per * cos_inclination ) alphaz = sin_arg_per * sin_inclination alpha = numpy . array ( [ alphax , alphay , alphaz ] ) betax = ( - cos_long_asc_nodes * sin_arg_per - sin_long_asc_nodes * cos_arg_per * cos_inclination ) betay = ( - sin_long_asc_nodes * sin_arg_per + cos_long_asc_nodes * cos_arg_per * cos_inclination ) betaz = cos_arg_per * sin_inclination beta = numpy . array ( [ betax , betay , betaz ] ) separation = ( semi_major_axis * ( 1.0 - eccentricity ** 2 ) / ( 1.0 + eccentricity * cos_true_anomaly ) ) position_vector = ( separation * cos_true_anomaly * alpha + separation * sin_true_anomaly * beta ) . T velocity_tilde = ( ( G * ( primary_mass + secondary_mass ) / ( semi_major_axis * ( 1.0 - eccentricity ** 2 ) ) ) ** 0.5 ) velocity_vector = ( - 1.0 * velocity_tilde * sin_true_anomaly * alpha + velocity_tilde * ( eccentricity + cos_true_anomaly ) * beta ) . T return position_vector , velocity_vector
Returns relative positions/velocities for secondaries orbiting primaries.
Returns relative positions / velocities for secondaries orbiting primaries .
"""
    Returns relative positions/velocities for secondaries orbiting primaries.
    If primary_mass is a scalar, assumes the same primary for all secondaries.
    """ # arrays need to be equal to number of secondaries, or have just one value # alpha is a unit vector directed along the line of node # beta is a unit vector perpendicular to alpha and the orbital angular # momentum vector # Relative position and velocity # Compute the relative separation # Common factor
primary_mass
secondary_mass
semi_major_axis
eccentricity
true_anomaly
inclination
longitude_of_the_ascending_node
argument_of_periapsis
G