Law: displacement of falling body (#845)

* Add convenience function for vectors

* Create law outline and test

* Refactor vector addition and subtraction

* Add law description

* Fix vector arithmetic functions

* Fix vector addition function usage

* Add note about Big O

symplyphysics/core/vectors/arithmetics.py

@@ -38,36 +38,53 @@ def equal_vectors(vector_left: Vector, vector_right: Vector) -> bool:
     return True
 
 
-# Sum of two vectors
-# Sum of two vectors can be seen as a diagonal of the parallelogram, where vectors are adjacent sides of this parallelogram.
-# To subtract vectors, multiply one of the vectors to -1 and add them.
+# Add any number of vectors
 #NOTE: adding two non-cartesian vectors is not a trivial task. We suggest to convert them to cartesian
 #      vectors, add them and convert back.
-def add_cartesian_vectors(vector_left: Vector, vector_right: Vector) -> Vector:
-    if vector_left.coordinate_system != vector_right.coordinate_system:
-        raise TypeError(
-            f"Different coordinate systems in vectors: {str(vector_left.coordinate_system)} vs {str(vector_right.coordinate_system)}"
-        )
-    if vector_left.coordinate_system.coord_system_type != CoordinateSystem.System.CARTESIAN:
-        coord_name_from = CoordinateSystem.system_to_transformation_name(
-            vector_left.coordinate_system.coord_system_type)
-        raise ValueError(
-            f"Addition is only supported for cartesian coordinates: got {coord_name_from}")
-    if vector_right.coordinate_system.coord_system_type != CoordinateSystem.System.CARTESIAN:
-        coord_name_from = CoordinateSystem.system_to_transformation_name(
-            vector_right.coordinate_system.coord_system_type)
-        raise ValueError(
-            f"Addition is only supported for cartesian coordinates: got {coord_name_from}")
-    (list_left_extended, list_right_extended) = _extend_two_vectors(vector_left, vector_right)
-    result = list(
-        map(lambda lr: sympify(lr[0] + lr[1]), zip(list_left_extended, list_right_extended)))
-    return Vector(result, vector_left.coordinate_system)
-
-
-def subtract_cartesian_vectors(vector_left: Vector, vector_right: Vector) -> Vector:
+def add_cartesian_vectors(*vectors: Vector) -> Vector:
+    def add_two_cartesian_vectors(vector_left: Vector, vector_right: Vector) -> Vector:
+        if vector_left.coordinate_system != vector_right.coordinate_system:
+            raise ValueError(
+                "Vectors must have the same coordinate system, "
+                f"got {vector_left.coordinate_system} and {vector_right.coordinate_system}"
+            )
+
+        left_type = vector_left.coordinate_system.coord_system_type
+        right_type = vector_right.coordinate_system.coord_system_type
+
+        if left_type != CoordinateSystem.System.CARTESIAN or right_type != CoordinateSystem.System.CARTESIAN:
+            left_name = CoordinateSystem.system_to_transformation_name(left_type)
+            right_name = CoordinateSystem.system_to_transformation_name(right_type)
+            raise ValueError(f"Expect two Cartesian vectors, got {left_name} and {right_name}")
+
+        (list_left_extended, list_right_extended) = _extend_two_vectors(vector_left, vector_right)
+        result = [
+            sympify(left_component + right_component)
+            for (left_component, right_component)
+            in zip(list_left_extended, list_right_extended)
+        ]
+
+        return Vector(result, vector_left.coordinate_system)
+
+    if len(vectors) < 1:
+        raise ValueError("Provide at least one vector")
+
+    result = vectors[0]
+    for vector in vectors[1:]:
+        result = add_two_cartesian_vectors(result, vector)
+    return result
+
+
+# Subtract vectors from the first given vector
+def subtract_cartesian_vectors(*vectors: Vector) -> Vector:
+    if len(vectors) < 2:
+        raise ValueError("Provide at least two vectors")
+
+    vector_minuend = vectors[0]
+    vector_subtrahend = add_cartesian_vectors(*vectors[1:])
     return add_cartesian_vectors(
-        vector_left,
-        scale_vector(-1, vector_right)
+        vector_minuend,
+        scale_vector(-1, vector_subtrahend),
     )
 
 

symplyphysics/laws/gravity/vector/falling_body_displacement.py

@@ -0,0 +1,103 @@
+from sympy import Expr, Rational
+from symplyphysics import (
+    angle_type,
+    units,
+    Quantity,
+    Vector,
+    QuantityVector,
+    validate_input,
+    validate_output,
+    scale_vector,
+    add_cartesian_vectors,
+    cross_cartesian_vectors,
+)
+
+# Description
+## Suppose a reference frame S' is fixed to a moving object A (e.g., Earth) and some body B moving freely
+## (i.e. the sum of external non-gravitational forces acting on it is zero). In the case of an inertial
+## frame of reference, the displacement of body B from the starting position would follow the usual rule
+## `s = v0 * t + (g / 2) * t**2`. But in the case of non-inertial frames, we have to take the Coriolis and
+## the centrifugal force into account as well, which results into the following series:
+
+# Conditions
+## - The sum `F` of all other, non-gravitational forces acting on body B is 0.
+## - `g` is independent of coordinates.
+
+# Notes
+## - Note that the series is truncated at the fifth term. More terms can be obtained by plugging the result
+##   into the equation of motion `a = g + cross(v, w)` and integrating it over time.
+
+# Law: s = v0 * t
+#        + t**2 * (g / 2 + cross(v0, w))
+#        + t**3 / 3 * (cross(g, w) + 2 * cross(cross(v0, w), w))
+#        + t**4 / 6 * cross(cross(g, w), w)
+#        + O(t**5)
+## t - time
+## s - vector of displacement of body B in S'
+## v0 - vector of initial velocity of body B
+## w - pseudovector of angular velocity of rotation of moving frame S'
+## g - vector of acceleration due to gravity (or free fall acceleration) of body B
+## cross(a, b) - cross product between vectors `a` and `b`
+## O - [big O](https://en.wikipedia.org/wiki/Big_O_notation) as per the formal mathematical definition
+##     in the limit of t approaching infinity
+
+
+def displacement_law(
+    time_: Expr,
+    initial_velocity_: Vector,
+    angular_velocity_: Vector,
+    acceleration_due_to_gravity_: Vector,
+) -> Vector:
+    initial_cross_angular = cross_cartesian_vectors(
+        initial_velocity_,
+        angular_velocity_,
+    )
+
+    acceleration_cross_angular = cross_cartesian_vectors(
+        acceleration_due_to_gravity_,
+        angular_velocity_,
+    )
+
+    return add_cartesian_vectors(
+        scale_vector(time_, initial_velocity_),
+        scale_vector(
+            time_**2,
+            add_cartesian_vectors(
+                scale_vector(Rational(1, 2), acceleration_due_to_gravity_),
+                initial_cross_angular,
+            ),
+        ),
+        scale_vector(
+            time_**3 / 3,
+            add_cartesian_vectors(
+                acceleration_cross_angular,
+                scale_vector(2, cross_cartesian_vectors(initial_cross_angular, angular_velocity_)),
+            ),
+        ),
+        scale_vector(
+            time_**4 / 6,
+            cross_cartesian_vectors(acceleration_cross_angular, angular_velocity_)
+        )
+    )
+
+
+@validate_input(
+    time_=units.time,
+    initial_velocity_=units.velocity,
+    angular_velocity_=angle_type / units.time,
+    acceleration_due_to_gravity_=units.acceleration,
+)
+@validate_output(units.length)
+def calculate_displacement(
+    time_: Quantity,
+    initial_velocity_: QuantityVector,
+    angular_velocity_: QuantityVector,
+    acceleration_due_to_gravity_: QuantityVector,
+) -> QuantityVector:
+    result_vector = displacement_law(
+        time_,
+        initial_velocity_.to_base_vector(),
+        angular_velocity_.to_base_vector(),
+        acceleration_due_to_gravity_.to_base_vector(),
+    )
+    return QuantityVector.from_base_vector(result_vector)

test/core/vectors/arithmetics_test.py

@@ -92,7 +92,7 @@ def test_rebased_add_vectors(test_args: Args) -> None:
 
 def test_invalid_add_vectors(test_args: Args) -> None:
     C1 = CoordinateSystem()
-    with raises(TypeError):
+    with raises(ValueError):
         add_cartesian_vectors(Vector([1, 2], test_args.C), Vector([1, 2], C1))
     # non-cartesian addition is not supported
     C2 = CoordinateSystem(CoordinateSystem.System.CYLINDRICAL)

test/gravity/vector/falling_body_displacement_test.py

@@ -0,0 +1,88 @@
+from collections import namedtuple
+from pytest import fixture, raises
+from symplyphysics import (
+    errors,
+    units,
+    Quantity,
+    QuantityVector,
+    assert_equal_vectors,
+)
+from symplyphysics.laws.gravity.vector import falling_body_displacement as law
+
+Args = namedtuple("Args", "t v0 w g")
+
+
+@fixture(name="test_args")
+def test_args_fixture() -> Args:
+    t = Quantity(5 * units.second)
+
+    v_unit = units.meter / units.second
+    v0 = QuantityVector([v_unit, v_unit, v_unit])
+
+    w = QuantityVector([0, -7.27e-5 * units.radian / units.second, 0])
+
+    g = QuantityVector([0, 0, -9.81 * units.meter / units.second**2])
+
+    return Args(t=t, v0=v0, w=w, g=g)
+
+
+def test_law(test_args: Args) -> None:
+    result = law.calculate_displacement(test_args.t, test_args.v0, test_args.w, test_args.g)
+    assert_equal_vectors(
+        result,
+        QuantityVector([4.97 * units.meter, 5.00 * units.meter, -118.0 * units.meter]),
+        tolerance=4e-3,
+    )
+
+
+def test_bad_time(test_args: Args) -> None:
+    tb = Quantity(units.coulomb)
+    with raises(errors.UnitsError):
+        law.calculate_displacement(tb, test_args.v0, test_args.w, test_args.g)
+    with raises(TypeError):
+        law.calculate_displacement(100, test_args.v0, test_args.w, test_args.g)
+
+
+def test_bad_linear_velocity(test_args: Args) -> None:
+    vb_vector = QuantityVector([units.coulomb])
+    with raises(errors.UnitsError):
+        law.calculate_displacement(test_args.t, vb_vector, test_args.w, test_args.g)
+
+    vb_scalar = Quantity(units.meter / units.second)
+    with raises(AttributeError):
+        law.calculate_displacement(test_args.t, vb_scalar, test_args.w, test_args.g)
+
+    with raises(TypeError):
+        law.calculate_displacement(test_args.t, 100, test_args.w, test_args.g)
+    with raises(TypeError):
+        law.calculate_displacement(test_args.t, [100], test_args.w, test_args.g)
+
+
+def test_bad_angular_velocity(test_args: Args) -> None:
+    wb_vector = QuantityVector([units.coulomb])
+    with raises(errors.UnitsError):
+        law.calculate_displacement(test_args.t, test_args.v0, wb_vector, test_args.g)
+
+    wb_scalar = Quantity(units.radian / units.second)
+    with raises(AttributeError):
+        law.calculate_displacement(test_args.t, test_args.v0, wb_scalar, test_args.g)
+
+    with raises(TypeError):
+        law.calculate_displacement(test_args.t, test_args.v0, 100, test_args.g)
+    with raises(TypeError):
+        law.calculate_displacement(test_args.t, test_args.v0, [100], test_args.g)
+
+
+def test_bad_acceleration(test_args: Args) -> None:
+    gb_vector = QuantityVector([units.coulomb])
+    with raises(errors.UnitsError):
+        law.calculate_displacement(test_args.t, test_args.v0, test_args.w, gb_vector)
+
+    gb_scalar = Quantity(units.meter / units.second**2)
+    with raises(AttributeError):
+        law.calculate_displacement(test_args.t, test_args.v0, test_args.w, gb_scalar)
+
+    with raises(TypeError):
+        law.calculate_displacement(test_args.t, test_args.v0, test_args.w, 100)
+    with raises(TypeError):
+        law.calculate_displacement(test_args.t, test_args.v0, test_args.w, [100])