Once you are at the python command line the first step is to import basic functionality from SymPy and the Mechanics module, otherwise you will only have basic python commands available to work with. We will import the symbols
function from SymPy core and with the * method
bring in all functionality from the mechanics package.
In [1]:
from sympy import symbols
from sympy.physics.mechanics import *
You can now see what functions and variables that are available to you with::
In [2]:
dir()
Out[2]:
This is a long list of available functions. Read about the python import statement to learn about better ways to import only what you need. One good explanation is http://effbot.org/zone/import-confusion.htm.
To get started working with vectors we will need to create a reference frame, as all vectors need to be defined with respect to a reference frame in the mechanics package. If you know the name of the command that you want to use simply use the builtin help function to bring up the documentation for the function. In our case we need to use the ReferenceFrame class.
In [3]:
help(ReferenceFrame)
The ReferenceFrame
class manages everything about rotations, angular velocities, and angular accelerations with respect to other reference frames. Now create an inertial reference frame called N for "Newtonian" as was described in the the docstring:
In [4]:
N = ReferenceFrame('N')
Keep in mind that N
is the variable name of which the reference frame named "N"
is stored. It is important to note that N
is an object and it has properties
and functions associated with it. To see a list of them type:
In [5]:
dir(N)
Out[5]:
Notice that three of the properties are x
, y
, and z
. These are the
orthonormal unit vectors associated with the reference frame and are the
building blocks for creating more general vectors. We can create a vector by simply
building a linear combination of the unit vectors:
In [6]:
v = 1 * N.x + 2 * N.y + 3 * N.z
Now a vector expressed in the N reference frame is stored in the variable v
.
We can print v
to the screen by typing:
In [7]:
print(v)
The vector v
can be manipulated as expected. You can multiply and divide them by scalars:
In [8]:
2 * v
Out[8]:
In [9]:
v / 3.0
Out[9]:
Note that three is expressed as 3.0
instead of 3
. The python language does
integer division by default. There are ways around this, but for now simply
remember to always declare numbers as floats (i.e. include a decimal).
You can also add and subtract vectors::
In [10]:
v + v
Out[10]:
In [11]:
w = 5 * N.x + 7 * N.y
In [12]:
v - w
Out[12]:
Vectors also have some useful properties:
In [13]:
dir(v)
Out[13]:
You can find the magnitude of a vector by typing:
In [14]:
v.magnitude()
Out[14]:
You can compute a unit vector in the direction of v
:
In [15]:
v.normalize()
Out[15]:
You can find the measure numbers and the reference frame the vector was defined in with::
In [16]:
v.args
Out[16]:
Dot and cross products of vectors can also be computed::
In [17]:
dot(v, w)
Out[17]:
In [18]:
cross(v, w)
Out[18]:
We've only used numbers as our measure numbers so far, but it is just as easy
to use symbols. We will introduce six symbols for our measure numbers with the
SymPy symbols
[help(symbols) for the documentation
] function:
In [19]:
a1, a2, a3 = symbols('a1 a2 a3')
b1, b2, b3 = symbols('b1 b2 b3')
And create two new vectors that are completely symbolic:
In [20]:
x = a1 * N.x + a2 * N.y + a3 * N.z
y = b1 * N.x + b2 * N.y + b3 * N.z
In [21]:
dot(x, y)
Out[21]:
In [22]:
z = cross(x, y)
z
Out[22]:
Numbers and symbols work together seamlessly:
In [23]:
dot(v, x)
Out[23]:
You can also differentiate a vector with respect to a variable in a reference frame:
In [24]:
z.diff(a1, N)
Out[24]:
Vectors don't have be defined with respect to just one reference frame. We can
create a new reference frame and orient it with respect to the N
frame that
has already been created. We will use the orient
method of the new frame to
do a simple rotation through alpha
about the N.x
axis:
In [25]:
A = ReferenceFrame('A')
alpha = symbols('alpha')
A.orient(N, 'Axis', [alpha, N.x])
Now the direction cosine matrix with of A
with respect to N
can be
computed:
In [26]:
A.dcm(N)
Out[26]:
Now that SymPy knows that A
and N
are oriented with respect to each other
we can express the vectors that we originally wrote in the A
frame:
In [27]:
v.express(A)
Out[27]:
In [28]:
z.express(A)
Out[28]:
In dynamics systems, at least some of the relative orientation of reference
frames and vectors are time varying. The mechanics module provides a way to
specify quantities as time varying. Let's define two variables beta
and d
as
variables which are functions of time:
In [29]:
beta, d = dynamicsymbols('beta d')
beta, d
Out[29]:
Now we can create a new reference frame that is oriented with respect to the A
frame by beta
and create a vector in that new frame that is a function of d
.
This time we will use the orientnew
method of the A
frame to create the new
reference frame B
:
In [30]:
B = A.orientnew('B', 'Axis', (beta, A.y))
In [31]:
vec = d * B.z
We can now compute the angular velocity of the reference frame B
with respect
to other reference frames:
In [32]:
B.ang_vel_in(N)
Out[32]:
This allows us to now differentiate the vector, vec
, with respect to time and
a reference frame::
In [33]:
vecdot = vec.dt(N)
vecdot
Out[33]:
In [34]:
vecdot.express(N)
Out[34]:
The dynamicsymbols
function also allows you to easily create the derivatives of time
varying variables and store them in a variable.
In [35]:
theta = dynamicsymbols('theta')
thetad = dynamicsymbols('theta', 1)
theta, thetad
Out[35]:
At this point we now have all the tools needed to setup the kinematics for a dynamic system. Let's start with the classic mass spring damper system under the influence of gravity. Go to the next notebook to follow along: https://github.com/PythonDynamics/pydy_examples/tree/master/mass_spring_damper.