Frames
Frame A
$$ {}^{N}{}{\vec{\omega}}^{A}{} = \dot{\theta} \hat{n}_z= \dot{\theta} \hat{a}_z $$
Velocity
if point $p$’s position is represented by $\vec{r}$ where
$$ \vec{r} = -l\hat{a}_y $$ then from the example in the derivatives module, we know that $$ \vec{v} = \frac{ {}^{N}{}{ d\vec{r}}{}{}}{dt}= l\dot{\theta}\hat{a}_x $$
and $$ \vec{a} = l\ddot{\theta}\hat{a}_x + l\dot{\theta}^2\hat{a}_y $$
Speed Variables
In this system our speed variable is $\dot{\theta}$
Forces
There are a number of forces acting within the system
- Two tension forces, one at point $p$ one at point $o$. These may be considered internal forces if modeled correctly
- One gravity force acting on point $p$
Partial Velocities
- partial velocity of particle
- $v_{p1} = v$
- partial velocity of f1, f2,f3
- $v_{f1} = 0$
- $v_{f2} = v$
- $v_{f3} = v$
Generalized effective force of a particle $$ {}^{N}{}{F}{Q}{u_r} = \frac{\partial{^N\vec{v}^Q}}{\partial u_r} \cdot \left(m^Q* {}^{N}{}{\vec{a}}^{Q}{}\right) $$
Kane’s Equation for a system of particles and masses for speed variable $u_r$, $m$ forces, $n$ moments, $o$ particles, and $p$ bodies $$ \sum_{i=1}^{m}{\vec{f}i} \cdot \frac{\partial\vec{v}{\vec{f}i}}{\partial u_r} + \sum{j=1}^{n}{\vec{\tau}j} \cdot \frac{\partial\vec{\omega}{\vec{\tau}j}}{\partial u_r} = \sum{k=1}^{o}\frac{\partial\vec{v}{k}}{ \partial u_r} \cdot m{k}\vec{a}{k} + \sum{l=1}^{p}\left( \frac{\partial\vec{v}{l}}{\partial u_r} \cdot m{l}\vec{a}{l} + \frac{\partial\vec{\omega}{l}}{\partial u_r} \cdot \left( \overrightarrow{I}_{l}\cdot\vec{\alpha}_l
- \vec{\omega}l \times \overrightarrow{I}{l}\cdot\vec{\omega}_l \right) \right) $$
Example 1
This example uses $\vec{r}=l\hat{a}_y$ as the position vector for describing particle p’s position
%matplotlib inline
# -*- coding: utf-8 -*-
"""
Written by Daniel M. Aukes
Email: danaukes<at>gmail.com
Please see LICENSE for full license.
"""
import pynamics
from pynamics.frame import Frame
from pynamics.variable_types import Differentiable,Constant,Variable
from pynamics.system import System
from pynamics.body import Body
from pynamics.dyadic import Dyadic
from pynamics.output import Output,PointsOutput
from pynamics.particle import Particle
import pynamics.integration
import sympy
import numpy
import matplotlib.pyplot as plt
plt.rcParams["animation.html"] = "html5"
plt.ion()
from math import pi
from matplotlib import animation, rc
from IPython.display import HTML
system = System()
pynamics.set_system(__name__,system)
lA = Constant(1,'lA',system)
mA = Constant(1,'mA',system)
g = Constant(9.81,'g',system)
tinitial = 0
tfinal = 5
fps = 30
t = numpy.r_[tinitial:tfinal:1/fps]
tol=1e-12
qA,qA_d,qA_dd = Differentiable('qA',system)
initialvalues = {}
initialvalues[qA]=0*pi/180
initialvalues[qA_d]=0*pi/180
ini = [initialvalues[item] for item in system.get_state_variables()]
N = Frame('N')
A = Frame('A')
system.set_newtonian(N)
A.rotate_fixed_axis_directed(N,[0,0,1],qA,system)
pO = 0*N.x
pNA=lA*A.x
ParticleA = Particle(pNA,mA,'ParticleA',system)
system.addforcegravity(-g*N.y)
f,ma = system.getdynamics()
2020-12-16 08:35:57,486 - pynamics.system - INFO - getting dynamic equations
f
[-g⋅lA⋅mA⋅cos(qA)]
ma
⎡ 2 ⎤
⎣lA ⋅mA⋅qA_dd⎦
func = system.state_space_post_invert(f,ma)
2020-12-16 08:35:57,571 - pynamics.system - INFO - solving a = f/m and creating function
2020-12-16 08:35:57,571 - pynamics.system - INFO - substituting constrained in Ma-f.
2020-12-16 08:35:57,671 - pynamics.system - INFO - done solving a = f/m and creating function
states=pynamics.integration.integrate(func,ini,t,rtol=tol,atol=tol,args=({'constants':system.constant_values},))
2020-12-16 08:35:57,972 - pynamics.integration - INFO - beginning integration
2020-12-16 08:35:57,972 - pynamics.system - INFO - integration at time 0000.00
2020-12-16 08:35:58,219 - pynamics.integration - INFO - finished integration
points = [pO,pNA]
points_output = PointsOutput(points,system)
points_output.calc(states)
points_output.animate(fps = fps,lw=2)
2020-12-16 08:35:58,241 - pynamics.output - INFO - calculating outputs
2020-12-16 08:35:58,250 - pynamics.output - INFO - done calculating outputs
del A,N,ParticleA,system, qA,qA_d,qA_dd,ini,states,t,tfinal,tinitial,tol
del pO,pNA
del lA,g,mA
del f,ma
del points
Example 2
This example uses $\vec{r}= x\hat{n}_x + y\hat{n}_y$ as the position vector for p. This requires an added constraint to maintain the length of r, as in $|\vec{r}|=l_A$
system = System()
pynamics.set_system(__name__,system)
lA = Constant(1,'lA',system)
mA = Constant(1,'mA',system)
g = Constant(9.81,'g',system)
tinitial = 0
tfinal = 5
fps = 30
t = numpy.r_[tinitial:tfinal:1/fps]
tol=1e-12
x,x_d,x_dd = Differentiable('x',system)
y,y_d,y_dd = Differentiable('y',system)
initialvalues = {}
initialvalues[x]=1
initialvalues[x_d]=0
initialvalues[y]=0
initialvalues[y_d]=0
ini = [initialvalues[item] for item in system.get_state_variables()]
N = Frame('N')
system.set_newtonian(N)
pO = 0*N.x
pNA=x*N.x+y*N.y
ParticleA = Particle(pNA,mA,'ParticleA',system)
system.addforcegravity(-g*N.y)
eq = []
eq.append(pNA.length() - lA)
eq_d=[(system.derivative(item)) for item in eq]
eq_dd=[(system.derivative(item)) for item in eq_d]
f,ma = system.getdynamics()
2020-12-16 08:35:58,873 - pynamics.system - INFO - getting dynamic equations
f
[0, -g⋅mA]
ma
[mA⋅x_dd, mA⋅y_dd]
func = system.state_space_post_invert(f,ma,eq_dd)
2020-12-16 08:35:58,958 - pynamics.system - INFO - solving a = f/m and creating function
2020-12-16 08:35:58,958 - pynamics.system - INFO - substituting constrained in Ma-f.
2020-12-16 08:35:59,021 - pynamics.system - INFO - done solving a = f/m and creating function
states=pynamics.integration.integrate(func,ini,t,rtol=tol,atol=tol,args=({'constants':system.constant_values},))
2020-12-16 08:35:59,043 - pynamics.integration - INFO - beginning integration
2020-12-16 08:35:59,043 - pynamics.system - INFO - integration at time 0000.00
2020-12-16 08:35:59,290 - pynamics.system - INFO - integration at time 0003.90
2020-12-16 08:35:59,359 - pynamics.integration - INFO - finished integration
points = [pO,pNA]
points_output2 = PointsOutput(points,system)
points_output2.calc(states)
points_output2.animate(fps = fps,lw=2)
2020-12-16 08:35:59,376 - pynamics.output - INFO - calculating outputs
2020-12-16 08:35:59,388 - pynamics.output - INFO - done calculating outputs
HTML(points_output.anim.to_html5_video())
HTML(points_output2.anim.to_html5_video())