How to use the poliastro.util.norm function in poliastro

To help you get started, we’ve selected a few poliastro examples, based on popular ways it is used in public projects.

Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.

github poliastro / poliastro / tests / tests_twobody / test_propagation.py View on Github external
def test_propagating_to_certain_nu_is_correct():
    # take an elliptic orbit
    a = 1.0 * u.AU
    ecc = 1.0 / 3.0 * u.one
    _a = 0.0 * u.rad
    nu = 10 * u.deg
    elliptic = Orbit.from_classical(Sun, a, ecc, _a, _a, _a, nu)

    elliptic_at_perihelion = elliptic.propagate_to_anomaly(0.0 * u.rad)
    r_per, _ = elliptic_at_perihelion.rv()

    elliptic_at_aphelion = elliptic.propagate_to_anomaly(np.pi * u.rad)
    r_ap, _ = elliptic_at_aphelion.rv()

    assert_quantity_allclose(norm(r_per), a * (1.0 - ecc))
    assert_quantity_allclose(norm(r_ap), a * (1.0 + ecc))

    # TODO: Test specific values
    assert elliptic_at_perihelion.epoch > elliptic.epoch
    assert elliptic_at_aphelion.epoch > elliptic.epoch

    # test 10 random true anomaly values
    # TODO: Rework this test
    for nu in np.random.uniform(low=-np.pi, high=np.pi, size=10):
        elliptic = elliptic.propagate_to_anomaly(nu * u.rad)
        r, _ = elliptic.rv()
        assert_quantity_allclose(norm(r), a * (1.0 - ecc ** 2) / (1 + ecc * np.cos(nu)))
github poliastro / poliastro / tests / tests_twobody / test_propagation.py View on Github external
def test_propagation_hyperbolic():
    # Data from Curtis, example 3.5
    r0 = [Earth.R.to(u.km).value + 300, 0, 0] * u.km
    v0 = [0, 15, 0] * u.km / u.s
    expected_r_norm = 163180 * u.km
    expected_v_norm = 10.51 * u.km / u.s

    ss0 = Orbit.from_vectors(Earth, r0, v0)
    tof = 14941 * u.s
    ss1 = ss0.propagate(tof)
    r, v = ss1.rv()

    assert_quantity_allclose(norm(r), expected_r_norm, rtol=1e-4)
    assert_quantity_allclose(norm(v), expected_v_norm, rtol=1e-3)
github poliastro / poliastro / src / poliastro / threebody / restricted.py View on Github external
), "Body 1 is not the main body: it has less mass than the 'secondary' body"

    # Define local frame of reference:
    # Center: main body, NOT the barycenter
    # x-axis: points to the secondary body
    ux = r2 - r1
    r12 = norm(ux)
    ux = ux / r12

    # y-axis: contained in the orbital plane, perpendicular to x-axis

    def cross(x, y):
        return np.cross(x, y) * x.unit * y.unit

    uy = cross(n, ux)
    uy = uy / norm(uy)

    # position in x-axis
    x1, x2, x3, x4, x5 = lagrange_points(r12, m1, m2)

    # position in y-axis
    # L1, L2, L3 are located in the x-axis, so y123 = 0

    # L4 and L5 are points in the plane of rotation which form an equilateral
    # triangle with the two masses (Battin)
    # sqrt(3)/2 = sin(60 deg)
    y4 = np.sqrt(3) / 2 * r12
    y5 = -y4

    # Convert L points coordinates (x,y) to original vectorial base [r1 r2]
    L1 = r1 + ux * x1
    L2 = r1 + ux * x2
github poliastro / poliastro / src / poliastro / plotting / tisserand.py View on Github external
num_contours: int
            Number of contour lines for flyby speed
        N: int
            Number of points for flyby angle

        Note
        ----
        The algorithm for generating Tisserand plots is the one depicted in
        "Preliminary Trajectory Design of a Mission to Enceladus" by David
        Falcato Fialho Palma, section 3.6

        """

        # Generate mean orbital elements Earth
        body_rv = get_mean_elements(body).to_vectors()
        R_body, V_body = norm(body_rv.r), norm(body_rv.v)

        # Generate non-dimensional velocity and alpha span
        vinf_array = np.linspace(vinf_span[0], vinf_span[-1], num_contours)
        alpha_array = np.linspace(alpha_lim[0], alpha_lim[-1], N)
        vinf_array /= V_body

        # Construct the mesh for any configuration
        V_INF, ALPHA = np.meshgrid(vinf_array, alpha_array)

        # Solving for non-dimensional a_sc and ecc_sc
        A_SC = 1 / np.abs(1 - V_INF ** 2 - 2 * V_INF * np.cos(ALPHA))
        ECC_SC = np.sqrt(1 - 1 / A_SC * ((3 - 1 / A_SC - V_INF ** 2) / (2)) ** 2)

        # Compute main Tisserand variables
        RR_P = A_SC * R_body * (1 - ECC_SC)
        RR_A = A_SC * R_body * (1 + ECC_SC)
github poliastro / poliastro / src / poliastro / plotting / porkchop.py View on Github external
ss_dpt = Orbit.from_vectors(attractor, rr_dpt_body, vv_dpt_body, epoch=t_launch)
    ss_arr = Orbit.from_vectors(attractor, rr_arr_body, vv_arr_body, epoch=t_arrival)

    # Define time of flight
    tof = ss_arr.epoch - ss_dpt.epoch

    if tof <= 0:
        return None, None, None, None, None

    try:
        # Lambert is now a Maneuver object
        man_lambert = Maneuver.lambert(ss_dpt, ss_arr)

        # Get norm delta velocities
        dv_dpt = norm(man_lambert.impulses[0][1])
        dv_arr = norm(man_lambert.impulses[1][1])

        # Compute all the output variables
        c3_launch = dv_dpt ** 2
        c3_arrival = dv_arr ** 2

        return (
            dv_dpt.to(u.km / u.s).value,
            dv_arr.to(u.km / u.s).value,
            c3_launch.to(u.km ** 2 / u.s ** 2).value,
            c3_arrival.to(u.km ** 2 / u.s ** 2).value,
            tof.jd,
        )

    except AssertionError:
        return None, None, None, None, None
github poliastro / poliastro / src / poliastro / threebody / flybys.py View on Github external
T_vec = np.cross(S_vec, c_vec) * u.one
    T_vec = T_vec / norm(T_vec)
    R_vec = np.cross(S_vec, T_vec) * u.one

    # This vector defines the B-Plane
    B_vec = b * (np.cos(theta) * T_vec + np.sin(theta) * R_vec)

    # We have to rotate the inbound hyperbolic excess velocity
    # an angle delta (turn angle) around a vector that is orthogonal to
    # the B-Plane and trajectory plane
    rot_v = np.cross(B_vec / b, S_vec) * u.one

    # And now we rotate the outbound hyperbolic excess velocity
    # u_vec = v_inf_1 / norm(v_inf) = S_vec
    v_vec = np.cross(rot_v, v_inf_1) * u.one
    v_vec = v_vec / norm(v_vec)

    v_inf_2 = v_inf * (np.cos(delta) * S_vec + np.sin(delta) * v_vec)

    v_spacecraft_out = v_inf_2 + v_body

    return v_spacecraft_out, delta
github poliastro / poliastro / src / poliastro / threebody / flybys.py View on Github external
v_inf = norm(v_inf_1)

    ecc = 1 + r_p * v_inf ** 2 / k  # Eccentricity of the entry hyperbola
    delta = 2 * np.arcsin(1 / ecc)  # Turn angle

    b = k / v_inf ** 2 * np.sqrt(ecc ** 2 - 1)  # Magnitude of the B vector

    # Now we compute the unit vectors in which to return the outbound hyperbolic excess velocity:
    # * S goes along the hyperbolic excess velocity and is perpendicular to the B-Plane,
    # * T goes along the B-Plane and is parallel to _some_ fundamental plane - in this case, the plane in which
    #   the velocities are computed
    # * R completes the orthonormal set
    S_vec = v_inf_1 / v_inf
    c_vec = np.array([0, 0, 1]) * u.one
    T_vec = np.cross(S_vec, c_vec) * u.one
    T_vec = T_vec / norm(T_vec)
    R_vec = np.cross(S_vec, T_vec) * u.one

    # This vector defines the B-Plane
    B_vec = b * (np.cos(theta) * T_vec + np.sin(theta) * R_vec)

    # We have to rotate the inbound hyperbolic excess velocity
    # an angle delta (turn angle) around a vector that is orthogonal to
    # the B-Plane and trajectory plane
    rot_v = np.cross(B_vec / b, S_vec) * u.one

    # And now we rotate the outbound hyperbolic excess velocity
    # u_vec = v_inf_1 / norm(v_inf) = S_vec
    v_vec = np.cross(rot_v, v_inf_1) * u.one
    v_vec = v_vec / norm(v_vec)

    v_inf_2 = v_inf * (np.cos(delta) * S_vec + np.sin(delta) * v_vec)
github poliastro / poliastro / src / poliastro / iod / izzo.py View on Github external
def _lambert(k, r1, r2, tof, M, numiter, rtol):
    # Check preconditions
    assert tof > 0
    assert k > 0

    # Chord
    c = r2 - r1
    c_norm, r1_norm, r2_norm = norm(c), norm(r1), norm(r2)

    # Semiperimeter
    s = (r1_norm + r2_norm + c_norm) * .5

    # Versors
    i_r1, i_r2 = r1 / r1_norm, r2 / r2_norm
    i_h = np.cross(i_r1, i_r2)
    i_h = i_h / norm(i_h)  # Fixed from paper

    # Geometry of the problem
    ll = np.sqrt(1 - c_norm / s)

    if i_h[2] < 0:
        ll = -ll
        i_h = -i_h
github poliastro / poliastro / src / poliastro / plotting.py View on Github external
def set_frame(self, p_vec, q_vec, w_vec):
        """Sets perifocal frame.

        Raises
        ------
        ValueError
            If the vectors are not a set of mutually orthogonal unit vectors.
        """
        if not np.allclose([norm(v) for v in (p_vec, q_vec, w_vec)], 1):
            raise ValueError("Vectors must be unit.")
        elif not np.allclose([p_vec.dot(q_vec), q_vec.dot(w_vec), w_vec.dot(p_vec)], 0):
            raise ValueError("Vectors must be mutually orthogonal.")
        else:
            self._frame = p_vec, q_vec, w_vec

        if self._orbits:
            self._redraw()