/

Source code for sympy.physics.mechanics.essential

__all__ = ['ReferenceFrame', 'Vector', 'Dyadic', 'dynamicsymbols',
'MechanicsStrPrinter', 'MechanicsPrettyPrinter',
'MechanicsLatexPrinter']

from sympy import (Matrix, Symbol, sin, cos, eye, trigsimp, diff, sqrt, sympify,
expand, zeros, Derivative, Function, symbols, Add,
solve, S)
from sympy.core import C
from sympy.core.function import UndefinedFunction
from sympy.printing.conventions import split_super_sub
from sympy.printing.latex import LatexPrinter
from sympy.printing.pretty.pretty import PrettyPrinter
from sympy.printing.pretty.stringpict import prettyForm, stringPict
from sympy.printing.str import StrPrinter
from sympy.utilities import group
from sympy.core.compatibility import reduce

[docs] def express(self, frame1, frame2=None): """Expresses this Dyadic in alternate frame(s) The first frame is the list side expression, the second frame is the right side; if Dyadic is in form A.x|B.y, you can express it in two different frames. If no second frame is given, the Dyadic is expressed in only one frame. Parameters ========== frame1 : ReferenceFrame The frame to express the left side of the Dyadic in frame2 : ReferenceFrame If provided, the frame to express the right side of the Dyadic in Examples ======== >>> from sympy.physics.mechanics import ReferenceFrame, outer, dynamicsymbols >>> N = ReferenceFrame('N') >>> q = dynamicsymbols('q') >>> B = N.orientnew('B', 'Axis', [q, N.z]) >>> d = outer(N.x, N.x) >>> d.express(B, N) cos(q)*(B.x|N.x) - sin(q)*(B.y|N.x) """ if frame2 == None: frame2 = frame1 _check_frame(frame1) _check_frame(frame2) ol = S(0) for i, v in enumerate(self.args): ol += v[0] * (v[1].express(frame1) | v[2].express(frame2)) return ol
[docs] def doit(self, **hints): """Calls .doit() on each term in the Dyadic""" return sum([Dyadic( [ (v[0].doit(**hints), v[1], v[2]) ]) for v in self.args])
[docs] def dt(self, frame): """Take the time derivative of this Dyadic in a frame. Parameters ========== frame : ReferenceFrame The frame to take the time derivative in Examples ======== >>> from sympy.physics.mechanics import ReferenceFrame, outer, dynamicsymbols >>> N = ReferenceFrame('N') >>> q = dynamicsymbols('q') >>> B = N.orientnew('B', 'Axis', [q, N.z]) >>> d = outer(N.x, N.x) >>> d.dt(B) - q'*(N.y|N.x) - q'*(N.x|N.y) """ _check_frame(frame) t = dynamicsymbols._t ol = S(0) for i, v in enumerate(self.args): ol += (v[0].diff(t) * (v[1] | v[2])) ol += (v[0] * (v[1].dt(frame) | v[2])) ol += (v[0] * (v[1] | v[2].dt(frame))) return ol
[docs] def simplify(self): """Simplify the elements in the Dyadic in-place.""" for i, v in enumerate(self.args): self.args[i] = (v[0].simplify(), v[1], v[2])
[docs] def subs(self, *args, **kwargs): """Substituion on the Dyadic. Examples ======== >>> from sympy.physics.mechanics import ReferenceFrame >>> from sympy import Symbol >>> N = ReferenceFrame('N') >>> s = Symbol('s') >>> a = s * (N.x|N.x) >>> a.subs({s: 2}) 2*(N.x|N.x) """ return sum([ Dyadic([(v[0].subs(*args, **kwargs), v[1], v[2])]) for v in self.args])
dot = __and__ cross = __xor__
[docs]class ReferenceFrame(object): """A reference frame in classical mechanics. ReferenceFrame is a class used to represent a reference frame in classical mechanics. It has a standard basis of three unit vectors in the frame's x, y, and z directions. It also can have a rotation relative to a parent frame; this rotation is defined by a direction cosine matrix relating this frame's basis vectors to the parent frame's basis vectors. It can also have an angular velocity vector, defined in another frame. """ def __init__(self, name, indices=None, latexs=None): """ReferenceFrame initialization method. A ReferenceFrame has a set of orthonormal basis vectors, along with orientations relative to other ReferenceFrames and angular velocities relative to other ReferenceFrames. Parameters ========== indices : list (of strings) If custom indices are desired for console, pretty, and LaTeX printing, supply three as a list. The basis vectors can then be accessed with the get_item method. latexs : list (of strings) If custom names are desired for LaTeX printing of each basis vector, supply the names here in a list. Examples ======== >>> from sympy.physics.mechanics import ReferenceFrame, mlatex >>> N = ReferenceFrame('N') >>> N.x N.x >>> O = ReferenceFrame('O', ('1', '2', '3')) >>> O.x O['1'] >>> O['1'] O['1'] >>> P = ReferenceFrame('P', latexs=('A1', 'A2', 'A3')) >>> mlatex(P.x) 'A1' """ if not isinstance(name, (str, unicode)): raise TypeError('Need to supply a valid name') # The if statements below are for custom printing of basis-vectors for # each frame. # First case, when custom indices are supplied if indices != None: if not isinstance(indices, (tuple, list)): raise TypeError('Supply the indices as a list') if len(indices) != 3: raise ValueError('Supply 3 indices') for i in indices: if not isinstance(i, (str, unicode)): raise TypeError('Indices must be strings') self.str_vecs = [(name + '[\'' + indices[0] + '\']'), (name + '[\'' + indices[1] + '\']'), (name + '[\'' + indices[2] + '\']')] self.pretty_vecs = [(u"\033[94m\033[1m" + name.lower() + u"_" + indices[0] + u"\033[0;0m\x1b[0;0m"), (u"\033[94m\033[1m" + name.lower() + u"_" + indices[1] + u"\033[0;0m\x1b[0;0m"), (u"\033[94m\033[1m" + name.lower() + u"_" + indices[2] + u"\033[0;0m\x1b[0;0m")] self.latex_vecs = [(r"\mathbf{\hat{%s}_{%s}}" % (name.lower(), indices[0])), (r"\mathbf{\hat{%s}_{%s}}" % (name.lower(), indices[1])), (r"\mathbf{\hat{%s}_{%s}}" % (name.lower(), indices[2]))] self.indices = indices # Second case, when no custom indices are supplied else: self.str_vecs = [(name + '.x'), (name + '.y'), (name + '.z')] self.pretty_vecs = [(u"\033[94m\033[1m" + name.lower() + u"_x\033[0;0m\x1b[0;0m"), (u"\033[94m\033[1m" + name.lower() + u"_y\033[0;0m\x1b[0;0m"), (u"\033[94m\033[1m" + name.lower() + u"_z\033[0;0m\x1b[0;0m")] self.latex_vecs = [(r"\mathbf{\hat{%s}_x}" % name.lower()), (r"\mathbf{\hat{%s}_y}" % name.lower()), (r"\mathbf{\hat{%s}_z}" % name.lower())] self.indices = ['x', 'y', 'z'] # Different step, for custom latex basis vectors if latexs != None: if not isinstance(latexs, (tuple, list)): raise TypeError('Supply the indices as a list') if len(latexs) != 3: raise ValueError('Supply 3 indices') for i in latexs: if not isinstance(i, (str, unicode)): raise TypeError('Latex entries must be strings') self.latex_vecs = latexs self.name = name self._dcm_dict = {} self._ang_vel_dict = {} self._ang_acc_dict = {} self._dlist = [self._dcm_dict, self._ang_vel_dict, self._ang_acc_dict] self._cur = 0 self._x = Vector([(Matrix([1, 0, 0]), self)]) self._y = Vector([(Matrix([0, 1, 0]), self)]) self._z = Vector([(Matrix([0, 0, 1]), self)]) def __getitem__(self, ind): """Returns basis vector for the provided index (index being an str)""" if not isinstance(ind, (str, unicode)): raise TypeError('Supply a valid str for the index') if self.indices[0] == ind: return self.x if self.indices[1] == ind: return self.y if self.indices[2] == ind: return self.z else: raise ValueError('Not a defined index') def __iter__(self): return iter([self.x, self.y, self.z]) def __str__(self): """Returns the name of the frame. """ return self.name __repr__ = __str__ def _dict_list(self, other, num): """Creates a list from self to other using _dcm_dict. """ outlist = [[self]] oldlist = [[]] while outlist != oldlist: oldlist = outlist[:] for i, v in enumerate(outlist): templist = v[-1]._dlist[num].keys() for i2, v2 in enumerate(templist): if not v.__contains__(v2): littletemplist = v + [v2] if not outlist.__contains__(littletemplist): outlist.append(littletemplist) for i, v in enumerate(oldlist): if v[-1] != other: outlist.remove(v) outlist.sort(key = len) if len(outlist) != 0: return outlist[0] raise ValueError('No Connecting Path found between ' + self.name + ' and ' + other.name) def _w_diff_dcm(self, otherframe): """Angular velocity from time differentiating the DCM. """ dcm2diff = self.dcm(otherframe) diffed = dcm2diff.diff(dynamicsymbols._t) angvelmat = diffed * dcm2diff.T w1 = trigsimp(expand(angvelmat[7]), recursive=True) w2 = trigsimp(expand(angvelmat[2]), recursive=True) w3 = trigsimp(expand(angvelmat[3]), recursive=True) return -Vector([(Matrix([w1, w2, w3]), self)])
[docs] def ang_acc_in(self, otherframe): """Returns the angular acceleration Vector of the ReferenceFrame. Effectively returns the Vector: ^N alpha ^B which represent the angular acceleration of B in N, where B is self, and N is otherframe. Parameters ========== otherframe : ReferenceFrame The ReferenceFrame which the angular acceleration is returned in. Examples ======== >>> from sympy.physics.mechanics import ReferenceFrame, Vector >>> N = ReferenceFrame('N') >>> A = ReferenceFrame('A') >>> V = 10 * N.x >>> A.set_ang_acc(N, V) >>> A.ang_acc_in(N) 10*N.x """ _check_frame(otherframe) if otherframe in self._ang_acc_dict: return self._ang_acc_dict[otherframe] else: return self.ang_vel_in(otherframe).dt(otherframe)
[docs] def ang_vel_in(self, otherframe): """Returns the angular velocity Vector of the ReferenceFrame. Effectively returns the Vector: ^N omega ^B which represent the angular velocity of B in N, where B is self, and N is otherframe. Parameters ========== otherframe : ReferenceFrame The ReferenceFrame which the angular velocity is returned in. Examples ======== >>> from sympy.physics.mechanics import ReferenceFrame, Vector >>> N = ReferenceFrame('N') >>> A = ReferenceFrame('A') >>> V = 10 * N.x >>> A.set_ang_vel(N, V) >>> A.ang_vel_in(N) 10*N.x """ _check_frame(otherframe) flist = self._dict_list(otherframe, 1) outvec = Vector([]) for i in range(len(flist) - 1): outvec += flist[i]._ang_vel_dict[flist[i + 1]] return outvec
[docs] def dcm(self, otherframe): """The direction cosine matrix between frames. This gives the DCM between this frame and the otherframe. The format is N.xyz = N.dcm(B) * B.xyz A SymPy Matrix is returned. Parameters ========== otherframe : ReferenceFrame The otherframe which the DCM is generated to. Examples ======== >>> from sympy.physics.mechanics import ReferenceFrame, Vector >>> from sympy import symbols >>> q1 = symbols('q1') >>> N = ReferenceFrame('N') >>> A = N.orientnew('A', 'Axis', [q1, N.x]) >>> N.dcm(A) [1, 0, 0] [0, cos(q1), -sin(q1)] [0, sin(q1), cos(q1)] """ _check_frame(otherframe) flist = self._dict_list(otherframe, 0) outdcm = eye(3) for i in range(len(flist) - 1): outdcm = outdcm * flist[i + 1]._dcm_dict[flist[i]] return outdcm
[docs] def orient(self, parent, rot_type, amounts, rot_order=''): """Defines the orientation of this frame relative to a parent frame. Supported orientation types are Body, Space, Quaternion, Axis. Examples show correct usage. Parameters ========== parent : ReferenceFrame The frame that this ReferenceFrame will have its orientation matrix defined in relation to. rot_type : str The type of orientation matrix that is being created. amounts : list OR value The quantities that the orientation matrix will be defined by. rot_order : str If applicable, the order of a series of rotations. Examples ======== >>> from sympy.physics.mechanics import ReferenceFrame, Vector >>> from sympy import symbols >>> q0, q1, q2, q3, q4 = symbols('q0 q1 q2 q3 q4') >>> N = ReferenceFrame('N') >>> B = ReferenceFrame('B') Now we have a choice of how to implement the orientation. First is Body. Body orientation takes this reference frame through three successive simple rotations. Acceptable rotation orders are of length 3, expressed in XYZ or 123, and cannot have a rotation about about an axis twice in a row. >>> B.orient(N, 'Body', [q1, q2, q3], '123') >>> B.orient(N, 'Body', [q1, q2, 0], 'ZXZ') >>> B.orient(N, 'Body', [0, 0, 0], 'XYX') Next is Space. Space is like Body, but the rotations are applied in the opposite order. >>> B.orient(N, 'Space', [q1, q2, q3], '312') Next is Quaternion. This orients the new ReferenceFrame with Quaternions, defined as a finite rotation about lambda, a unit vector, by some amount theta. This orientation is described by four parameters: q0 = cos(theta/2) q1 = lambda_x sin(theta/2) q2 = lambda_y sin(theta/2) q3 = lambda_z sin(theta/2) Quaternion does not take in a rotation order. >>> B.orient(N, 'Quaternion', [q0, q1, q2, q3]) Last is Axis. This is a rotation about an arbitrary, non-time-varying axis by some angle. The axis is supplied as a Vector. This is how simple rotations are defined. >>> B.orient(N, 'Axis', [q1, N.x + 2 * N.y]) """ _check_frame(parent) amounts = list(amounts) for i, v in enumerate(amounts): if not isinstance(v, Vector): amounts[i] = sympify(v) def _rot(axis, angle): """DCM for simple axis 1,2,or 3 rotations. """ if axis == 1: return Matrix([[1, 0, 0], [0, cos(angle), -sin(angle)], [0, sin(angle), cos(angle)]]) elif axis == 2: return Matrix([[cos(angle), 0, sin(angle)], [0, 1, 0], [-sin(angle), 0, cos(angle)]]) elif axis == 3: return Matrix([[cos(angle), -sin(angle), 0], [sin(angle), cos(angle), 0], [0, 0, 1]]) approved_orders = ('123', '231', '312', '132', '213', '321', '121', '131', '212', '232', '313', '323', '') rot_order = str(rot_order).upper() # Now we need to make sure XYZ = 123 rot_type = rot_type.upper() rot_order = [i.replace('X', '1') for i in rot_order] rot_order = [i.replace('Y', '2') for i in rot_order] rot_order = [i.replace('Z', '3') for i in rot_order] rot_order = ''.join(rot_order) if not rot_order in approved_orders: raise TypeError('The supplied order is not an approved type') parent_orient = [] if rot_type == 'AXIS': if not rot_order == '': raise TypeError('Axis orientation takes no rotation order') if not (isinstance(amounts, (list, tuple)) & (len(amounts) == 2)): raise TypeError('Amounts are a list or tuple of length 2') theta = amounts[0] axis = amounts[1] axis = _check_vector(axis) if not axis.dt(parent) == 0: raise ValueError('Axis cannot be time-varying') axis = axis.express(parent).normalize() axis = axis.args[0][0] parent_orient = ((eye(3) - axis * axis.T) * cos(theta) + Matrix([[0, -axis[2], axis[1]],[axis[2], 0, -axis[0]], [-axis[1], axis[0], 0]]) * sin(theta) + axis * axis.T) elif rot_type == 'QUATERNION': if not rot_order == '': raise TypeError('Quaternion orientation takes no rotation order') if not (isinstance(amounts, (list, tuple)) & (len(amounts) == 4)): raise TypeError('Amounts are a list or tuple of length 4') q0, q1, q2, q3 = amounts parent_orient = (Matrix([[q0 ** 2 + q1 ** 2 - q2 ** 2 - q3 ** 2, 2 * (q1 * q2 - q0 * q3), 2 * (q0 * q2 + q1 * q3)], [2 * (q1 * q2 + q0 * q3), q0 ** 2 - q1 ** 2 + q2 **2 - q3 ** 2, 2 * (q2 * q3 - q0 * q1)], [2 * (q1 * q3 - q0 * q2), 2 * (q0 * q1 + q2 * q3), q0 ** 2 - q1 ** 2 - q2 ** 2 + q3 ** 2]])) elif rot_type == 'BODY': if not (len(amounts) == 3 & len(rot_order) == 3): raise TypeError('Body orientation takes 3 values & 3 orders') a1 = int(rot_order[0]) a2 = int(rot_order[1]) a3 = int(rot_order[2]) parent_orient = (_rot(a1, amounts[0]) * _rot(a2, amounts[1]) * _rot(a3, amounts[2])) elif rot_type == 'SPACE': if not (len(amounts) == 3 & len(rot_order) == 3): raise TypeError('Space orientation takes 3 values & 3 orders') a1 = int(rot_order[0]) a2 = int(rot_order[1]) a3 = int(rot_order[2]) parent_orient = (_rot(a3, amounts[2]) * _rot(a2, amounts[1]) * _rot(a1, amounts[0])) else: raise NotImplementedError('That is not an implemented rotation') self._dcm_dict.update({parent: parent_orient}) parent._dcm_dict.update({self: parent_orient.T}) if rot_type == 'QUATERNION': t = dynamicsymbols._t q0, q1, q2, q3 = amounts q0d = diff(q0, t) q1d = diff(q1, t) q2d = diff(q2, t) q3d = diff(q3, t) w1 = 2 * (q1d * q0 + q2d * q3 - q3d * q2 - q0d * q1) w2 = 2 * (q2d * q0 + q3d * q1 - q1d * q3 - q0d * q2) w3 = 2 * (q3d * q0 + q1d * q2 - q2d * q1 - q0d * q3) wvec = Vector([(Matrix([w1, w2, w3]), self)]) elif rot_type == 'AXIS': thetad = (amounts[0]).diff(dynamicsymbols._t) wvec = thetad * amounts[1].express(parent).normalize() else: try: from sympy.polys.polyerrors import CoercionFailed from sympy.physics.mechanics.functions import kinematic_equations q1, q2, q3 = amounts u1, u2, u3 = dynamicsymbols('u1, u2, u3') templist = kinematic_equations([u1, u2, u3], [q1, q2, q3], rot_type, rot_order) templist = [expand(i) for i in templist] td = solve(templist, [u1, u2, u3]) u1 = expand(td[u1]) u2 = expand(td[u2]) u3 = expand(td[u3]) wvec = u1 * self.x + u2 * self.y + u3 * self.z except (CoercionFailed, AssertionError): wvec = self._w_diff_dcm(parent) self._ang_vel_dict.update({parent: wvec}) parent._ang_vel_dict.update({self: -wvec})
[docs] def orientnew(self, newname, rot_type, amounts, rot_order='', indices=None, latexs=None): """Creates a new ReferenceFrame oriented with respect to this Frame. See ReferenceFrame.orient() for acceptable rotation types, amounts, and orders. Parent is going to be self. Parameters ========== newname : str The name for the new ReferenceFrame rot_type : str The type of orientation matrix that is being created. amounts : list OR value The quantities that the orientation matrix will be defined by. rot_order : str If applicable, the order of a series of rotations. Examples ======== >>> from sympy.physics.mechanics import ReferenceFrame, Vector >>> from sympy import symbols >>> q1 = symbols('q1') >>> N = ReferenceFrame('N') >>> A = N.orientnew('A', 'Axis', [q1, N.x]) .orient() documentation:\n ======================== """ newframe = ReferenceFrame(newname, indices, latexs) newframe.orient(self, rot_type, amounts, rot_order) return newframe
orientnew.__doc__ += orient.__doc__
[docs] def set_ang_acc(self, otherframe, value): """Define the angular acceleration Vector in a ReferenceFrame. Defines the angular acceleration of this ReferenceFrame, in another. Angular acceleration can be defined with respect to multiple different ReferenceFrames. Care must be taken to not create loops which are inconsistent. Parameters ========== otherframe : ReferenceFrame A ReferenceFrame to define the angular acceleration in value : Vector The Vector representing angular acceleration Examples ======== >>> from sympy.physics.mechanics import ReferenceFrame, Vector >>> N = ReferenceFrame('N') >>> A = ReferenceFrame('A') >>> V = 10 * N.x >>> A.set_ang_acc(N, V) >>> A.ang_acc_in(N) 10*N.x """ value = _check_vector(value) _check_frame(otherframe) self._ang_acc_dict.update({otherframe: value}) otherframe._ang_acc_dict.update({self: -value})
[docs] def set_ang_vel(self, otherframe, value): """Define the angular velocity vector in a ReferenceFrame. Defines the angular velocity of this ReferenceFrame, in another. Angular velocity can be defined with respect to multiple different ReferenceFrames. Care must be taken to not create loops which are inconsistent. Parameters ========== otherframe : ReferenceFrame A ReferenceFrame to define the angular velocity in value : Vector The Vector representing angular velocity Examples ======== >>> from sympy.physics.mechanics import ReferenceFrame, Vector >>> N = ReferenceFrame('N') >>> A = ReferenceFrame('A') >>> V = 10 * N.x >>> A.set_ang_vel(N, V) >>> A.ang_vel_in(N) 10*N.x """ value = _check_vector(value) _check_frame(otherframe) self._ang_vel_dict.update({otherframe: value}) otherframe._ang_vel_dict.update({self: -value})
@property
[docs] def x(self): """The basis Vector for the ReferenceFrame, in the x direction. """ return self._x
@property
[docs] def y(self): """The basis Vector for the ReferenceFrame, in the y direction. """ return self._y
@property
[docs] def z(self): """The basis Vector for the ReferenceFrame, in the z direction. """ return self._z
[docs] def dot(self, other): return self & other
dot.__doc__ = __and__.__doc__
[docs] def cross(self, other): return self ^ other
cross.__doc__ = __xor__.__doc__
[docs] def outer(self, other): return self | other
outer.__doc__ = __or__.__doc__
[docs] def diff(self, wrt, otherframe): """Takes the partial derivative, with respect to a value, in a frame. Returns a Vector. Parameters ========== wrt : Symbol What the partial derivative is taken with respect to. otherframe : ReferenceFrame The ReferenceFrame that the partial derivative is taken in. Examples ======== >>> from sympy.physics.mechanics import ReferenceFrame, Vector, dynamicsymbols >>> from sympy import Symbol >>> Vector.simp = True >>> t = Symbol('t') >>> q1 = dynamicsymbols('q1') >>> N = ReferenceFrame('N') >>> A = N.orientnew('A', 'Axis', [q1, N.y]) >>> A.x.diff(t, N) - q1'*A.z """ wrt = sympify(wrt) _check_frame(otherframe) outvec = S(0) for i,v in enumerate(self.args): if v[1] == otherframe: outvec += Vector([(v[0].diff(wrt), otherframe)]) 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
[docs] def doit(self, **hints): """Calls .doit() on each term in the Vector""" ov = S(0) for i, v in enumerate(self.args): ov += Vector([(v[0].applyfunc(lambda x: x.doit(**hints)) , v[1])]) return ov
[docs] def dt(self, otherframe): """Returns the time derivative of the Vector in a ReferenceFrame. Returns a Vector which is the time derivative of the self Vector, taken in frame otherframe. Parameters ========== otherframe : ReferenceFrame The ReferenceFrame that the partial derivative is taken in. Examples ======== >>> from sympy.physics.mechanics import ReferenceFrame, Vector, dynamicsymbols >>> from sympy import Symbol >>> q1 = Symbol('q1') >>> u1 = dynamicsymbols('u1') >>> N = ReferenceFrame('N') >>> A = N.orientnew('A', 'Axis', [q1, N.x]) >>> v = u1 * N.x >>> A.set_ang_vel(N, 10*A.x) >>> A.x.dt(N) == 0 True >>> v.dt(N) u1'*N.x """ outvec = S(0) _check_frame(otherframe) for i, v in enumerate(self.args): if v[1] == otherframe: outvec += Vector([(v[0].diff(dynamicsymbols._t), otherframe)]) else: outvec += (Vector([v]).dt(v[1]) + (v[1].ang_vel_in(otherframe) ^ Vector([v]))) return outvec
[docs] def express(self, otherframe): """Returns a vector, expressed in the other frame. A new Vector is returned, equalivalent to this Vector, but its components are all defined in only the otherframe. Parameters ========== otherframe : ReferenceFrame The frame for this Vector to be described in Examples ======== >>> from sympy.physics.mechanics import ReferenceFrame, Vector, dynamicsymbols >>> q1 = dynamicsymbols('q1') >>> N = ReferenceFrame('N') >>> A = N.orientnew('A', 'Axis', [q1, N.y]) >>> A.x.express(N) cos(q1)*N.x - sin(q1)*N.z """ _check_frame(otherframe) outvec = Vector(self.args + []) for i, v in enumerate(self.args): if v[1] != otherframe: temp = otherframe.dcm(v[1]) * v[0] for i2, v2 in enumerate(temp): if Vector.simp == True: temp[i2] = trigsimp(v2, recursive=True) else: temp[i2] = v2 outvec += Vector([(temp, otherframe)]) outvec -= Vector([v]) return outvec
[docs] def simplify(self): """Simplify the elements in the Vector in place. """ for i in self.args: i[0].simplify()
[docs] def subs(self, *args, **kwargs): """Substituion on the Vector. Examples ======== >>> from sympy.physics.mechanics import ReferenceFrame >>> from sympy import Symbol >>> N = ReferenceFrame('N') >>> s = Symbol('s') >>> a = N.x * s >>> a.subs({s: 2}) 2*N.x """ ov = S(0) for i, v in enumerate(self.args): ov += Vector([(v[0].subs(*args, **kwargs), v[1])]) return ov
[docs] def magnitude(self): """Returns the magnitude (Euclidean norm) of self.""" return sqrt(self & self)
[docs] def normalize(self): """Returns a Vector of magnitude 1, codirectional with self.""" return Vector(self.args + []) / self.magnitude()
class MechanicsStrPrinter(StrPrinter): """String Printer for mechanics. """ def _print_Derivative(self, e): t = dynamicsymbols._t if (bool(sum([i == t for i in e.variables])) & isinstance(type(e.args[0]), UndefinedFunction)): ol = str(e.args[0].func) for i, v in enumerate(e.variables): ol += dynamicsymbols._str return ol else: return StrPrinter().doprint(e) def _print_Function(self, e): t = dynamicsymbols._t if isinstance(type(e), UndefinedFunction): return StrPrinter().doprint(e).replace("(%s)"%t, '') return e.func.__name__ + "(%s)"%self.stringify(e.args, ", ") class MechanicsLatexPrinter(LatexPrinter): """Latex Printer for mechanics. """ def _print_Function(self, expr, exp=None): func = expr.func.__name__ t = dynamicsymbols._t if hasattr(self, '_print_' + func): return getattr(self, '_print_' + func)(expr, exp) elif isinstance(type(expr), UndefinedFunction) and (expr.args == (t,)): name, sup, sub = split_super_sub(func) if len(sup) != 0: sup = r"^{%s}" % "".join(sup) else: sup = r"" if len(sub) != 0: sub = r"_{%s}" % "".join(sub) else: sub = r"" return r"%s" % (name + sup + sub) else: args = [ str(self._print(arg)) for arg in expr.args ] # How inverse trig functions should be displayed, formats are: # abbreviated: asin, full: arcsin, power: sin^-1 inv_trig_style = self._settings['inv_trig_style'] # If we are dealing with a power-style inverse trig function inv_trig_power_case = False # If it is applicable to fold the argument brackets can_fold_brackets = self._settings['fold_func_brackets'] and \ len(args) == 1 and \ not self._needs_function_brackets(expr.args[0]) inv_trig_table = ["asin", "acos", "atan", "acot"] # If the function is an inverse trig function, handle the style if func in inv_trig_table: if inv_trig_style == "abbreviated": func = func elif inv_trig_style == "full": func = "arc" + func[1:] elif inv_trig_style == "power": func = func[1:] inv_trig_power_case = True # Can never fold brackets if we're raised to a power if exp is not None: can_fold_brackets = False if inv_trig_power_case: name = r"\operatorname{%s}^{-1}" % func elif exp is not None: name = r"\operatorname{%s}^{%s}" % (func, exp) else: name = r"\operatorname{%s}" % func if can_fold_brackets: name += r"%s" else: name += r"\left(%s\right)" if inv_trig_power_case and exp is not None: name += r"^{%s}" % exp return name % ",".join(args) def _print_Derivative(self, der_expr): # make sure it is an the right form der_expr = der_expr.doit() if not isinstance(der_expr, Derivative): return self.doprint(der_expr) # check if expr is a dynamicsymbol from sympy.core.function import AppliedUndef t = dynamicsymbols._t expr = der_expr.expr red = expr.atoms(AppliedUndef) syms = der_expr.variables test1 = not all([True for i in red if i.atoms() == set([t])]) test2 = not all([(t == i) for i in syms]) if test1 or test2: return LatexPrinter().doprint(der_expr) # done checking dots = len(syms) base = self._print_Function(expr) base_split = base.split('_', 1) base = base_split[0] if dots == 1: base = r"\dot{%s}" % base elif dots == 2: base = r"\ddot{%s}" % base elif dots == 3: base = r"\dddot{%s}" % base if len(base_split) is not 1: base += '_' + base_split[1] return base class MechanicsPrettyPrinter(PrettyPrinter): """Pretty Printer for mechanics. """ def _print_Derivative(self, deriv): # XXX use U('PARTIAL DIFFERENTIAL') here ? t = dynamicsymbols._t dots = 0 can_break = True syms = list(reversed(deriv.variables)) x = None while len(syms) > 0: if syms[-1] == t: syms.pop() dots += 1 else: break f = prettyForm(binding=prettyForm.FUNC, *self._print(deriv.expr)) if not (isinstance(type(deriv.expr), UndefinedFunction) and (deriv.expr.args == (t,))): dots = 0 can_break = False f = prettyForm(binding=prettyForm.FUNC, *self._print(deriv.expr).parens()) if dots == 0: dots = u"" elif dots == 1: dots = u"\u0307" elif dots == 2: dots = u"\u0308" elif dots == 3: dots = u"\u20db" elif dots == 4: dots = u"\u20dc" uni_subs = [u"\u2080", u"\u2081", u"\u2082", u"\u2083", u"\u2084", u"\u2085", u"\u2086", u"\u2087", u"\u2088", u"\u2089", u"\u208a", u"\u208b", u"\u208c", u"\u208d", u"\u208e", u"\u208f", u"\u2090", u"\u2091", u"\u2092", u"\u2093", u"\u2094", u"\u2095", u"\u2096", u"\u2097", u"\u2098", u"\u2099", u"\u209a", u"\u209b", u"\u209c", u"\u209d", u"\u209e", u"\u209f"] fpic = f.__dict__['picture'] funi = f.__dict__['unicode'] ind = len(funi) val = "" for i in uni_subs: cur_ind = funi.find(i) if (cur_ind != -1) and (cur_ind < ind): ind = cur_ind val = i if ind == len(funi): funi += dots else: funi = funi.replace(val, dots + val) if f.__dict__['picture'] == [f.__dict__['unicode']]: fpic = [funi] f.__dict__['picture'] = fpic f.__dict__['unicode'] = funi if (len(syms)) == 0 and can_break: return f for sym, num in group(syms, multiple=False): s = self._print(sym) ds = prettyForm(*s.left('d')) if num > 1: ds = ds**prettyForm(str(num)) if x is None: x = ds else: x = prettyForm(*x.right(' ')) x = prettyForm(*x.right(ds)) pform = prettyForm('d') if len(syms) > 1: pform = pform**prettyForm(str(len(syms))) pform = prettyForm(*pform.below(stringPict.LINE, x)) pform.baseline = pform.baseline + 1 pform = prettyForm(*stringPict.next(pform, f)) return pform def _print_Function(self, e): t = dynamicsymbols._t # XXX works only for applied functions func = e.func args = e.args func_name = func.__name__ prettyFunc = self._print(C.Symbol(func_name)) prettyArgs = prettyForm(*self._print_seq(args).parens()) # If this function is an Undefined function of t, it is probably a # dynamic symbol, so we'll skip the (t). The rest of the code is # identical to the normal PrettyPrinter code if isinstance(func, UndefinedFunction) and (args == (t,)): pform = prettyForm(binding=prettyForm.FUNC, *stringPict.next(prettyFunc)) else: pform = prettyForm(binding=prettyForm.FUNC, *stringPict.next(prettyFunc, prettyArgs)) # store pform parts so it can be reassembled e.g. when powered pform.prettyFunc = prettyFunc pform.prettyArgs = prettyArgs return pform def _check_dyadic(other): if not isinstance(other, Dyadic): other = sympify(other) if other != S(0): raise TypeError('A Dyadic must be supplied') else: other = Dyadic([]) return other def _check_frame(other): if not isinstance(other, ReferenceFrame): raise TypeError('A ReferenceFrame must be supplied') def _check_vector(other): if not isinstance(other, Vector): other = sympify(other) if other != S(0): raise TypeError('A Vector must be supplied') else: other = Vector([]) return other
[docs]def dynamicsymbols(names, level=0): """Uses symbols and Function for functions of time. Creates a SymPy UndefinedFunction, which is then initialized as a function of a variable, the default being Symbol('t'). Parameters ========== names : str Names of the dynamic symbols you want to create; works the same way as inputs to symbols level : int Level of differentiation of the returned function; d/dt once of t, twice of t, etc. Examples ======= >>> from sympy.physics.mechanics import dynamicsymbols >>> from sympy import diff, Symbol >>> q1 = dynamicsymbols('q1') >>> q1 q1(t) >>> diff(q1, Symbol('t')) Derivative(q1(t), t) """ esses = symbols(names, cls=Function) t = dynamicsymbols._t if hasattr(esses, '__iter__'): esses = [reduce(diff, [t]*level, e(t)) for e in esses] return esses else: return reduce(diff, [t]*level, esses(t))
dynamicsymbols._t = Symbol('t') dynamicsymbols._str = '\''