Repository: sympy/sympy

Issue:
partial_velocity should not require a rexpression in the provided reference frame
I think that we may have a fundamental design flaw with some things about physics.vector. To illustrate, let's say that we have defined a vector that describes the linear velocity of point P in the reference frame N in terms of the generalized speeds: u1, u2, and u3.

```python
N_v_P = u1 * A.x + u2 * A.y + u3 * A.z
```

Note that we have this velocity vector is expressed in a frame A which has an arbitrary relationship with the frame N.

If point P happens to be an essential point to form the equations of motion with Kane's method, then we will need to take the partial velocity of the vector with respect to the generalized speeds in an intertial reference frame. If N happens to be that inertial reference frame, then we should be able to do:

```python
N_v_P.diff(u1, N)
```
and this should return:

```python
A.x
```
The partial velocity of N_v_P in N with respect to u1 is simply A.x. Of course, if A is oriented with respect to N in some way, then A.x can be expressed in N's unit vectors, but it isn't necessary to do that to compute the partial velocity.

The following code shows that this does not work:

```
moorepants@garuda:pydy(error-messages)$ ipython
Python 3.5.1 |Continuum Analytics, Inc.| (default, Dec  7 2015, 11:16:01) 
Type "copyright", "credits" or "license" for more information.

IPython 4.1.2 -- An enhanced Interactive Python.
?         -> Introduction and overview of IPython's features.
%quickref -> Quick reference.
help      -> Python's own help system.
object?   -> Details about 'object', use 'object??' for extra details.

In [1]: import sympy.physics.mechanics as me

In [2]: u1, u2, u3 = me.dynamicsymbols('u1, u2, u3')

In [3]: N = me.ReferenceFrame('N')

In [4]: A = me.ReferenceFrame('A')

In [5]: N_v_P = u1 * A.x + u2 * A.y + u3 * A.z

In [7]: N_v_P.diff(u1, N)
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-7-7a11f17c841e> in <module>()
----> 1 N_v_P.diff(u1, N)

/home/moorepants/miniconda3/lib/python3.5/site-packages/sympy/physics/vector/vector.py in diff(self, wrt, otherframe)
    519                 outvec += Vector([(v[0].diff(wrt), otherframe)])
    520             else:
--> 521                 if otherframe.dcm(v[1]).diff(wrt) == zeros(3, 3):
    522                     d = v[0].diff(wrt)
    523                     outvec += Vector([(d, v[1])])

/home/moorepants/miniconda3/lib/python3.5/site-packages/sympy/physics/vector/frame.py in dcm(self, otherframe)
    405         if otherframe in self._dcm_cache:
    406             return self._dcm_cache[otherframe]
--> 407         flist = self._dict_list(otherframe, 0)
    408         outdcm = eye(3)
    409         for i in range(len(flist) - 1):

/home/moorepants/miniconda3/lib/python3.5/site-packages/sympy/physics/vector/frame.py in _dict_list(self, other, num)
    250             return outlist[0]
    251         raise ValueError('No Connecting Path found between ' + self.name +
--> 252                          ' and ' + other.name)
    253 
    254     def _w_diff_dcm(self, otherframe):

ValueError: No Connecting Path found between N and A

In [9]: me.partial_velocity([N_v_P], [u1, u2, u3], N)
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-9-56e053551919> in <module>()
----> 1 me.partial_velocity([N_v_P], [u1, u2, u3], N)

/home/moorepants/miniconda3/lib/python3.5/site-packages/sympy/physics/vector/functions.py in partial_velocity(vel_list, u_list, frame)
    564         pvlist = []
    565         for j in u_list:
--> 566             vel = i.diff(j, frame)
    567             pvlist += [vel]
    568         list_of_pvlists += [pvlist]

/home/moorepants/miniconda3/lib/python3.5/site-packages/sympy/physics/vector/vector.py in diff(self, wrt, otherframe)
    519                 outvec += Vector([(v[0].diff(wrt), otherframe)])
    520             else:
--> 521                 if otherframe.dcm(v[1]).diff(wrt) == zeros(3, 3):
    522                     d = v[0].diff(wrt)
    523                     outvec += Vector([(d, v[1])])

/home/moorepants/miniconda3/lib/python3.5/site-packages/sympy/physics/vector/frame.py in dcm(self, otherframe)
    405         if otherframe in self._dcm_cache:
    406             return self._dcm_cache[otherframe]
--> 407         flist = self._dict_list(otherframe, 0)
    408         outdcm = eye(3)
    409         for i in range(len(flist) - 1):

/home/moorepants/miniconda3/lib/python3.5/site-packages/sympy/physics/vector/frame.py in _dict_list(self, other, num)
    250             return outlist[0]
    251         raise ValueError('No Connecting Path found between ' + self.name +
--> 252                          ' and ' + other.name)
    253 
    254     def _w_diff_dcm(self, otherframe):

ValueError: No Connecting Path found between N and A
```

The problem is that the Vector object has no awareness of that it is a velocity of some point in some reference frame. Because of this, it doesn't really make sense to ask it to differentiate itself in the reference frame N.

But our vector differentiation in frames does work. It works because we *require* that there be a per-specified orientation between reference frames A and N. If this orientation is specified the correct answer is returned:

```
In [10]: q1, q2, q3 = me.dynamicsymbols('q1, q2, q3')

In [11]: A.orient(N, 'Body', (q1, q2, q3), 'XYZ')

In [12]: N_v_P.diff(u1, N)
Out[12]: A.x
```

This requires re-expression of the vector in the N frame so that the measure numbers in the N frame of N_v_P can be differentiated in that frame. But we know that re-expression is not necessary to get the correct result.

This comes up when you have ignorable coordinates. For example, if you want to study a satellite that is moving in inertial reference frame N, you can often ignore the 6 coordinates associated with the primary degrees of freedom up front. Currently, our software does not allow us to easily do this because of the above mentioned limitation.

Some thoughts:

- Maybe we need to define differentiation of the position/orientatoins, velocity, and accelerations at only on the Point and RefernceFrame objects because these have the full information about a given velocity vector, i.e. those objects know what frame the vector is with respect to. The vector by itself is not aware of what frame it is wrt.
- Maybe we need a "null" orientation option that lets you connect frames A and N but not specify the generalized coordinates and direction cosine matrix that relates the two. Then the vector.diff() function can deal with this properly.
- The fact that Vectors can differentiate themselves without know what frame they are defined wrt, seems mathematically incorrect.
- It is important to note that re-expression of a vector in another frame is different that knowing what frame a velocity vector is associated with.
- Do we need a VelocityVector object that knows more info like the frame the velocity is wrt or the associated points?
- I think Autolev properly tracks this with their syntax, for example in `N_v_P> = u1 * A1> + u2 * A2> + u3 * A3>` both the left and right hand sides have semantic meaning. This defines a velocity vector in N of P with that vector expression. All of that info must be stored in a data structure behind the scenes. They must have different types of vectors.

@sympy/mechanics What do you all think about this? I feel like we have an odd hole here that may take more than a hack to properly fix.

tl;dr Differentiation of vectors wrt to generalized speeds should not require that orientations be set between the frames that are associated with the unit vectors in the velocity vector expression and the frame that the velocity is defined it.

Discussion:

Patch:
diff --git a/sympy/physics/vector/functions.py b/sympy/physics/vector/functions.py
index a5cdabfc25..ba5a963c32 100644
--- a/sympy/physics/vector/functions.py
+++ b/sympy/physics/vector/functions.py
@@ -517,12 +517,10 @@ def _process_vector_differential(vectdiff, condition, \
         return (acc, vel, kwargs['position'])
 
 
-def partial_velocity(vel_list, u_list, frame):
-    """Returns a list of partial velocities.
-
-    For a list of velocity or angular velocity vectors the partial derivatives
-    with respect to the supplied generalized speeds are computed, in the
-    specified ReferenceFrame.
+def partial_velocity(vel_vecs, gen_speeds, frame):
+    """Returns a list of partial velocities with respect to the provided
+    generalized speeds in the given reference frame for each of the supplied
+    velocity vectors.
 
     The output is a list of lists. The outer list has a number of elements
     equal to the number of supplied velocity vectors. The inner lists are, for
@@ -532,12 +530,13 @@ def partial_velocity(vel_list, u_list, frame):
     Parameters
     ==========
 
-    vel_list : list
-        List of velocities of Point's and angular velocities of ReferenceFrame's
-    u_list : list
-        List of independent generalized speeds.
+    vel_vecs : iterable
+        An iterable of velocity vectors (angular or linear).
+    gen_speeds : iterable
+        An iterable of generalized speeds.
     frame : ReferenceFrame
-        The ReferenceFrame the partial derivatives are going to be taken in.
+        The reference frame that the partial derivatives are going to be taken
+        in.
 
     Examples
     ========
@@ -549,24 +548,27 @@ def partial_velocity(vel_list, u_list, frame):
     >>> N = ReferenceFrame('N')
     >>> P = Point('P')
     >>> P.set_vel(N, u * N.x)
-    >>> vel_list = [P.vel(N)]
-    >>> u_list = [u]
-    >>> partial_velocity(vel_list, u_list, N)
+    >>> vel_vecs = [P.vel(N)]
+    >>> gen_speeds = [u]
+    >>> partial_velocity(vel_vecs, gen_speeds, N)
     [[N.x]]
 
     """
-    if not iterable(vel_list):
-        raise TypeError('Provide velocities in an iterable')
-    if not iterable(u_list):
-        raise TypeError('Provide speeds in an iterable')
-    list_of_pvlists = []
-    for i in vel_list:
-        pvlist = []
-        for j in u_list:
-            vel = i.diff(j, frame)
-            pvlist += [vel]
-        list_of_pvlists += [pvlist]
-    return list_of_pvlists
+
+    if not iterable(vel_vecs):
+        raise TypeError('Velocity vectors must be contained in an iterable.')
+
+    if not iterable(gen_speeds):
+        raise TypeError('Generalized speeds must be contained in an iterable')
+
+    vec_partials = []
+    for vec in vel_vecs:
+        partials = []
+        for speed in gen_speeds:
+            partials.append(vec.diff(speed, frame, var_in_dcm=False))
+        vec_partials.append(partials)
+
+    return vec_partials
 
 
 def dynamicsymbols(names, level=0):
diff --git a/sympy/physics/vector/vector.py b/sympy/physics/vector/vector.py
index e7577df17b..4fbef9f85f 100644
--- a/sympy/physics/vector/vector.py
+++ b/sympy/physics/vector/vector.py
@@ -485,24 +485,30 @@ def outer(self, other):
         return self | other
     outer.__doc__ = __or__.__doc__
 
-    def diff(self, wrt, otherframe):
-        """Takes the partial derivative, with respect to a value, in a frame.
-
-        Returns a Vector.
+    def diff(self, var, frame, var_in_dcm=True):
+        """Returns the partial derivative of the vector with respect to a
+        variable in the provided reference frame.
 
         Parameters
         ==========
-
-        wrt : Symbol
+        var : Symbol
             What the partial derivative is taken with respect to.
-        otherframe : ReferenceFrame
-            The ReferenceFrame that the partial derivative is taken in.
+        frame : ReferenceFrame
+            The reference frame that the partial derivative is taken in.
+        var_in_dcm : boolean
+            If true, the differentiation algorithm assumes that the variable
+            may be present in any of the direction cosine matrices that relate
+            the frame to the frames of any component of the vector. But if it
+            is known that the variable is not present in the direction cosine
+            matrices, false can be set to skip full reexpression in the desired
+            frame.
 
         Examples
         ========
 
-        >>> from sympy.physics.vector import ReferenceFrame, Vector, dynamicsymbols
         >>> from sympy import Symbol
+        >>> from sympy.physics.vector import dynamicsymbols, ReferenceFrame
+        >>> from sympy.physics.vector import Vector
         >>> Vector.simp = True
         >>> t = Symbol('t')
         >>> q1 = dynamicsymbols('q1')
@@ -510,24 +516,39 @@ def diff(self, wrt, otherframe):
         >>> A = N.orientnew('A', 'Axis', [q1, N.y])
         >>> A.x.diff(t, N)
         - q1'*A.z
+        >>> B = ReferenceFrame('B')
+        >>> u1, u2 = dynamicsymbols('u1, u2')
+        >>> v = u1 * A.x + u2 * B.y
+        >>> v.diff(u2, N, var_in_dcm=False)
+        B.y
 
         """
 
         from sympy.physics.vector.frame import _check_frame
-        wrt = sympify(wrt)
-        _check_frame(otherframe)
-        outvec = Vector(0)
-        for i, v in enumerate(self.args):
-            if v[1] == otherframe:
-                outvec += Vector([(v[0].diff(wrt), otherframe)])
+
+        var = sympify(var)
+        _check_frame(frame)
+
+        partial = Vector(0)
+
+        for vector_component in self.args:
+            measure_number = vector_component[0]
+            component_frame = vector_component[1]
+            if component_frame == frame:
+                partial += Vector([(measure_number.diff(var), frame)])
             else:
-                if otherframe.dcm(v[1]).diff(wrt) == zeros(3, 3):
-                    d = v[0].diff(wrt)
-                    outvec += Vector([(d, v[1])])
-                else:
-                    d = (Vector([v]).express(otherframe)).args[0][0].diff(wrt)
-                    outvec += Vector([(d, otherframe)]).express(v[1])
-        return outvec
+                # If the direction cosine matrix relating the component frame
+                # with the derivative frame does not contain the variable.
+                if not var_in_dcm or (frame.dcm(component_frame).diff(var) ==
+                                      zeros(3, 3)):
+                    partial += Vector([(measure_number.diff(var),
+                                        component_frame)])
+                else:  # else express in the frame
+                    reexp_vec_comp = Vector([vector_component]).express(frame)
+                    deriv = reexp_vec_comp.args[0][0].diff(var)
+                    partial += Vector([(deriv, frame)]).express(component_frame)
+
+        return partial
 
     def express(self, otherframe, variables=False):
         """


Test patch:
diff --git a/sympy/physics/vector/tests/test_functions.py b/sympy/physics/vector/tests/test_functions.py
index 9a30180ac0..f6e561b5c0 100644
--- a/sympy/physics/vector/tests/test_functions.py
+++ b/sympy/physics/vector/tests/test_functions.py
@@ -1,11 +1,10 @@
 from sympy import S, Integral, sin, cos, pi, sqrt, symbols
-from sympy.physics.vector import (Dyadic, Point, ReferenceFrame, \
-                                  Vector)
-from sympy.physics.vector import (cross, dot, express, \
-                                  time_derivative, kinematic_equations, \
-                                  outer, partial_velocity, \
-                                  get_motion_params)
-from sympy.physics.vector.functions import dynamicsymbols
+from sympy.physics.vector import Dyadic, Point, ReferenceFrame, Vector
+from sympy.physics.vector.functions import (cross, dot, express,
+                                            time_derivative,
+                                            kinematic_equations, outer,
+                                            partial_velocity,
+                                            get_motion_params, dynamicsymbols)
 from sympy.utilities.pytest import raises
 
 Vector.simp = True
@@ -465,3 +464,10 @@ def test_partial_velocity():
             [[- r*L.y, r*L.x, 0, L.x, cos(q2)*L.y - sin(q2)*L.z],
             [0, 0, 0, L.x, cos(q2)*L.y - sin(q2)*L.z],
             [L.x, L.y, L.z, 0, 0]])
+
+    # Make sure that partial velocities can be computed regardless if the
+    # orientation between frames is defined or not.
+    A = ReferenceFrame('A')
+    B = ReferenceFrame('B')
+    v = u4 * A.x + u5 * B.y
+    assert partial_velocity((v, ), (u4, u5), A) == [[A.x, B.y]]
diff --git a/sympy/physics/vector/tests/test_vector.py b/sympy/physics/vector/tests/test_vector.py
index 602872d9ad..5079c8fac4 100644
--- a/sympy/physics/vector/tests/test_vector.py
+++ b/sympy/physics/vector/tests/test_vector.py
@@ -1,7 +1,7 @@
 from sympy import symbols, pi, sin, cos, ImmutableMatrix as Matrix
-from sympy.physics.vector import ReferenceFrame, Vector, \
-     dynamicsymbols, dot
+from sympy.physics.vector import ReferenceFrame, Vector, dynamicsymbols, dot
 from sympy.abc import x, y, z
+from sympy.utilities.pytest import raises
 
 
 Vector.simp = True
@@ -54,6 +54,7 @@ def test_Vector():
     assert v1.separate() == {A: v1}
     assert v5.separate() == {A: x*A.x + y*A.y, B: z*B.z}
 
+
 def test_Vector_diffs():
     q1, q2, q3, q4 = dynamicsymbols('q1 q2 q3 q4')
     q1d, q2d, q3d, q4d = dynamicsymbols('q1 q2 q3 q4', 1)
@@ -121,6 +122,25 @@ def test_Vector_diffs():
     assert v4.diff(q3d, B) == B.x + q3 * N.x + N.y
 
 
+def test_vector_var_in_dcm():
+
+    N = ReferenceFrame('N')
+    A = ReferenceFrame('A')
+    B = ReferenceFrame('B')
+    u1, u2, u3, u4 = dynamicsymbols('u1 u2 u3 u4')
+
+    v = u1 * u2 * A.x + u3 * N.y + u4**2 * N.z
+
+    assert v.diff(u1, N, var_in_dcm=False) == u2 * A.x
+    assert v.diff(u1, A, var_in_dcm=False) == u2 * A.x
+    assert v.diff(u3, N, var_in_dcm=False) == N.y
+    assert v.diff(u3, A, var_in_dcm=False) == N.y
+    assert v.diff(u3, B, var_in_dcm=False) == N.y
+    assert v.diff(u4, N, var_in_dcm=False) == 2 * u4 * N.z
+
+    raises(ValueError, lambda: v.diff(u1, N))
+
+
 def test_vector_simplify():
     x, y, z, k, n, m, w, f, s, A = symbols('x, y, z, k, n, m, w, f, s, A')
     N = ReferenceFrame('N')