Yes, what you have there is correct (B2 = [M]^-1 * KM.linearize()[1])
You didn't substitute in your constants on the second matrix you
calculate (B2.subs({g:9.8,m:1,l:1})). Once you do this, B1 and B2 will
not evaluate as equal, but if you do (B1 - B2.subs).simplify(), you
get all zeros.-Gilbert On Jan 23, 8:11 pm, Gustavo <[email protected]> wrote: > Thanks for the reply. I did mean for the torque to be due to an > actuator. > > Regarding the second point -- I didn't follow entirely. > > I understand that the resulting A and B matrices should be evaluated > at the point of linearization. In my case, then, they will have no > unbound symbols. I do want an approximation for A and B in qudots = > A*qu + B*tau. So then I should use > > B2 = [M]^-1 * KM.linearize()[1] > > B1 and B2 resemble each other but they're still not identical, so I > think I might still be misunderstanding. > > updated gist:https://gist.github.com/1665634 > > On Jan 23, 8:14 pm, Gilbert gede <[email protected]> wrote: > > > > > > > > > Two problems here: > > First: you've added a torque term which applies to ReferenceFrame B, > > tau*A.z. There are 2 potential issues with this. One, if possible, > > you should describe the torque in the same ReferenceFrame in which the > > angular velocity of the "torqued" frame is described in. Here, > > B.ang_vel_in(N) is described in the N ReferenceFrame, u2 * N.z, so it > > would best to describe the torque in the N ReferenceFrame (tau * N.z). > > In this particular example, it shouldn't make much of a difference, > > but in general you want to have as few intermediate frames as possible > > between vectors you are working with (or rather, as simple a direction > > cosine matrix relating the two frame's basis vectors as possible). > > The second issue with the torque you've added is that you only applied > > it to ReferenceFrame B. If you meant for this to be some "externally" > > applied torque, than this would be appropriate. If instead you meant > > for there to be something like a spring/actuator between the two > > pendulum links, you need to apply an equal and opposite torque to > > ReferenceFrame A, -tau * N.z. > > > Second: the linearization process only linearizes the forcing vector. > > If our equations are in the form [M]qudot = f (M is mass matrix, qudot > > is derivative of state vector, f is forcing vector), what the > > linearization functions returns is f_linearized as it would fit in > > this equation: [M]qudot = [f_lin_A]qu + [f_lin_B] tau . Linearization > > of the mass matrix should not be necessary; in fact it is checked that > > it is safe to linearize (no time-varying forces, masses, lengths, etc. > > show up in it). So what you would then have is qudot = [M]^-1 > > [f_lin_a]qu + [M]^-1 [f_lin_B] tau . Note that M, f_lin_A, and > > f_lin_B should be evaluated at the point of linearization. > > I now realize there are two issues with this. One, I did not write > > this docstring very clearly, so I will go and rewrite it to make this > > more clear. Two, it should return a vector associated with the B > > matrix, as it is not clear if you had say, two unique torques, what is > > associated with each column of the B matrix. I'll go and open a PR for > > these changes now. > > > One extra note: Also, I'm going to change it so that Particles and > > RigidBodies have an associated name, and take all relevant parameters > > on initialization (very unclear errors are given if you don't have > > everything). The pydy.org examples will be changed accordingly. > > > -Gilbert > > > On Jan 23, 1:39 pm, Gustavo <[email protected]> wrote: > > > > As far as I can tell, I have correctly added a torque term at the "P" > > > joint. > > > The linearization produces the matrices of expected dimension. I > > > haven't checked that the "A" matrix is correct, but it looks like the > > > "B" matrix is not correct. > > > > Here's a gist that computes "B" using the linearize() method (stored > > > in B2), and also in the way that I expect B to be computed (stored in > > > B1)https://gist.github.com/1665634 > > > > On Jan 19, 3:07 pm, Gustavo <[email protected]> wrote: > > > > > I found this example [1], which is very helpful. All I need is a > > > > torque term at the "P" joint. > > > > Then when I linearize the equations of motion, I hope to get a 4-by-4 > > > > "A" matrix and a 4-by-1 "B" matrix. > > > > > [1]http://pydy.org/index.php?title=Double_Pendulum#Integration_with_Scipy > > > > > On Jan 18, 3:41 pm, Gustavo <[email protected]> wrote: > > > > > > I'm looking to use the physics.mechanics module to get the equations > > > > > of motions for adoublependulumactuated at the center joint, and get > > > > > the linearized dynamics about certain points. I know that this can be > > > > > done using the Lagrangian, but I'd like to try to use the mechanics > > > > > module. > > > > > > I see the examples here [1] but I am having trouble applying that to > > > > > thedoublependulum. I have only an introductory classical mechanics > > > > > background, so I'm hoping that the mechanics module will let me > > > > > assemble a dynamical system, impose some constraints, and then > > > > > generate the equations of motion. > > > > > > Thanks, > > > > > Gustavo > > > > > > [1]http://docs.sympy.org/dev/modules/physics/mechanics/examples.html -- You received this message because you are subscribed to the Google Groups "sympy" group. To post to this group, send email to [email protected]. To unsubscribe from this group, send email to [email protected]. For more options, visit this group at http://groups.google.com/group/sympy?hl=en.
