I had a sign error. The code works now. Btw, sympy.physics.mechanics is likely overkill for high school physics. You can use it to solve these kinds of problem, but the overhead is high.
Jason moorepants.info +01 530-601-9791 On Mon, Feb 16, 2015 at 10:54 AM, Jason Moore <moorepa...@gmail.com> wrote: > Here is one way to do it: > > https://gist.github.com/moorepants/858503aa180df60a7829 > > But solve is returning an empty list. I'm not sure why that isn't working. > > > Jason > moorepants.info > +01 530-601-9791 > > On Mon, Feb 16, 2015 at 10:28 AM, Ondřej Čertík <ondrej.cer...@gmail.com> > wrote: > >> Hi Petr, >> >> On Mon, Feb 16, 2015 at 2:56 AM, Petr Baudis <pa...@ucw.cz> wrote: >> > Hello! >> > >> > This is going to be a rather long post, so I'll first explain what I'm >> > trying to do - then list all the essentially unsuccessful ways I tried >> > so far. >> > >> > Basically, my long-term goal is to build a database of high-school level >> > physics equations (dynamics, energy, electromagnetism, waves, ...) and >> > then be able to answer queries like "what is A given X, Y, Z?" with >> > formulas (or even numerical results). I need to make this completely >> > automatic, without any human input on which equations are relevant or >> > how to massage them to apply them to a problem. >> >> I think this might work to some extent, but I suspect some manual >> tweaking might be necessary. >> >> > >> > However, as an initial toy problem, I'm simply considering: >> > >> > Let us release a mass object at height $x_0$. Given gravity >> > uniform acceleration $g$, what will its landing velocity be? >> > >> > Basically, I'm looking for a way to supply SymPy some basic kinematic >> > equations and perform a query such that I get something like >> > >> > sym.Eq(v, sym.sqrt(2 * g * x_0)) >> > >> > on output. >> > >> > >> > I have not been very successful so far, though! (Using SymPy-0.7.6.) >> > Either (i) I'm doing various trival novice mistakes, (ii) there are >> some bugs >> > or easy-to-add missing features in SymPy preventing me to do this, >> (iii) what >> > I want to do is fundamentally hard problem, or (iv) SymPy is in too >> early >> > development stage and there is a more appropriate symbolic computation >> system >> > to use for this right now. I'd appreciate any advice! >> > >> > >> > Here's what I tried so far: >> > >> > (A) Most naive and least flexible approach - define functions s(t) and >> v(t), >> > build some simple equations around them and then ask for definition >> > of v(t) given that s(t) = x_0: >> > >> > import sympy as sym >> > import sympy.physics.mechanics as me >> > from sympy import init_printing >> > >> > x_0, t, g = sym.symbols('x_0 t g') >> > ua = sym.symbols('ua', positive=True) # acceleration is >> function of time, uniform acceleration is not >> > s, v = sym.symbols('s v', cls=sym.Function) # functions of t >> > initial = {v(0): 0, s(t): x_0} >> > >> > velocity_ds_eq = sym.Eq(v(t), s(t) / t) >> > unifaccel_dv_eq = sym.Eq(ua, (v(t) - v(0))/2/t) >> > eqs = [eq.subs(initial) for eq in [velocity_ds_eq, >> unifaccel_dv_eq]] >> > print(sym.solve(eqs, v(t), implicit=True)) >> > # -> [] >> >> I think you have to specify two unknowns to solve for, if you have a >> system of two equations, e.g. change v(t) -> [v(t), t], as you did in >> your second try: >> >> In [15]: print solve(eqs, [v(t), t]) >> [(-sqrt(2)*sqrt(ua)*sqrt(x_0), -sqrt(2)*sqrt(x_0)/(2*sqrt(ua))), >> (sqrt(2)*sqrt(ua)*sqrt(x_0), sqrt(2)*sqrt(x_0)/(2*sqrt(ua)))] >> >> Then it works. >> >> Why do you use "implicit=True"? It seems it doesn't work with it. >> >> > >> > Here's the kick - if I apply three rather arbitrary changes to the >> > last code segment, I *do* get a solution I'm after (thanks to >> > @vramana for helping me out with this!): >> > >> > velocity_ds_eq = sym.Eq(v(t) * t, s(t)) >> > unifaccel_dv_eq = sym.Eq(ua * t, (v(t) - v(0))/2) >> > eqs = [eq.subs(initial) for eq in [velocity_ds_eq, >> unifaccel_dv_eq]] >> > print(sym.solve(eqs, [v(t), t], implicit=True)) >> > # -> [(-sqrt(2)*ua*sqrt(x_0/ua), -sqrt(2)*sqrt(x_0/ua)/2), >> (sqrt(2)*ua*sqrt(x_0/ua), sqrt(2)*sqrt(x_0/ua)/2)] >> > >> > However, I don't understand why is it required to do either of these >> > changes - that is, solving for [v(t),t] (hmm, maybe I do have an >> idea >> > but not sure) and especially converting /t to *t in both of the >> > equations (doing that only in the first one will yield >> > {v(t): x_0/t, t: x_0/(2*t*ua)} - I guess I could deal with that). >> > The latter two changes give me an impression of a brittle approach >> > that might break randomly in more complex cases than the toy >> problem. >> >> I don't think multiplying by "t" has anything to do with the result, >> see my comment above. >> >> After we find a solution that you like, I can look at the other >> problems you mentioned. >> >> Ondrej >> >> > >> > (Of course I would prefer to represent velocity and acceleration as >> > vectors etc., but I got stuck much sooner than that!) >> > >> > >> > (B) Inspired by "Taming math and physics using SymPy", representing >> > velocity as integrated acceleration and position as integrated >> > velocity. I originally aimed to make both functions of time: >> > >> > import sympy as sym >> > import sympy.physics.mechanics as me >> > from sympy import init_printing >> > >> > t, a, v_0, x_0 = sym.symbols('t a v_0 x_0') >> > v, x = sym.symbols('v x', cls=sym.Function) # functions of (t) >> > # initial = {x_0: 100, v_0: 0, a: -9.8} >> > initial = {v_0: 0} >> > >> > eqv = sym.Eq(v(t), v_0 + sym.integrate(a, (t, 0,t))) >> > eqx = sym.Eq(x(t), x_0 + sym.integrate(v(t), (t, 0,t))) >> > eql = sym.Eq(x(t), 0) >> > eqs = [eq.subs(initial) for eq in eqx,eqv,eql] >> > print(sym.solve(eqs, [v(t), t], implicit=True)) >> > # -> [] >> > >> > Is my expectation that the sym.solve() call should yield the >> > solution I'm after correct? >> > >> > >> > (C) More literal to "Taming math and physics using SymPy", representing >> > vi and xi by simple expressions instead of functions: >> > >> > import sympy as sym >> > import sympy.physics.mechanics as me >> > from sympy import init_printing >> > >> > t, a, v, v_0, x_0 = sym.symbols('t a v v_0 x_0') >> > initial = {v_0: 0} >> > >> > vi = v_0 + sym.integrate(a, (t, 0,t)) >> > xi = x_0 + sym.integrate(vi, (t, 0,t)) >> > xeq = sym.Eq(xi, 0) >> > veq = sym.Eq(vi, v) >> > print(sym.solve([xeq, veq], [v,t], implicit=True)) >> > # -> {v: (-a*(a*t**2 + 2*x_0)/2 + v_0**2)/v_0, t: -(a*t**2/2 + >> x_0)/v_0} >> > >> > Basically, I'm not sure how to phrase my query - especially how to >> > obtain a solution for v which is independent of t. I thought that >> is >> > what the exclude= parameter of sym.solve() is for, but using it has >> > no effect. >> > >> > This works, but xeq is completely arbitrary invention, which is why >> > such a solution is of little use to me (it's not clear to me how to >> > come up with such equations in a generic way): >> > >> > xeq = sym.Eq((v*v), ((v*v).expand() - 2*a*x).simplify()) >> > print(sym.solve(xeq, v)[1].subs(initial)) >> > # -> sqrt(2)*sqrt(-a*x_0) >> > >> > >> > (D) Starting with position function x(t) and expressing velocity and >> > acceleration as time derivations of that. This was actually my >> > first approach (when I still hoped things would be easy ;). >> > >> > I tried to leverage sympy.physics.mechanics for this. I spent >> > quite a lot of time going through various examples and docs but >> > I think I still didn't fully get it. I can define the point and >> > references, describe velocity and acceleration, and solve the >> > differential equation for x: >> > >> > import sympy as sym >> > import sympy.physics.mechanics as me >> > from sympy import init_printing >> > >> > N = me.ReferenceFrame('N') # world reference frame >> > o = me.Point('o') # reference point (origin) for any world >> positions >> > o.set_vel(N, 0) >> > x_0, g = sym.symbols('x_0 g') >> > >> > r = me.Point('r') >> > x = me.dynamicsymbols('x') >> > r.set_pos(o, x * N.y) >> > >> > t = sym.symbols('t') >> > r_v = sym.Derivative(x, t) >> > r.set_vel(N, r_v * N.y) >> > >> > r_aeq = me.dot(g * N.y - r.acc(N), N.y) >> > r_xeq = sym.dsolve(r_aeq, x, ics={x.subs(t, 0): x_0, >> x.diff(t).subs(t, 0): 0}) >> > # -> x(t) == C1 + C2*t + g*t**2/2 >> > # Huh, ics parameter did not work: >> > r_xeq = r_xeq.replace('C1', x_0).replace('C2', 0) >> > # -> x(t) == x_0 + g*t**2/2 >> > >> > # x(t) == 0 on landing >> > # (I would've hoped that [r_xeq, x] would work instead of this:) >> > r_landing_eq = r_xeq.rhs >> > >> > # This step is a little arbitrary for my taste already: >> > r_teq = sym.solve(r_landing_eq, t) >> > # -> [-sqrt(2)*sqrt(-x_0/g), sqrt(2)*sqrt(-x_0/g)] >> > r_teq = r_teq[1] >> > >> > But now I'm stuck - I need to re-express r_v using r_aeq and r_teq. >> > My few naive attempts mostly using (d)solve failed. How to proceed? >> > >> > >> > I'd be glad for any assistance or advice. I hope someone finished >> > reading up to here... ;-) >> > >> > -- >> > Petr Baudis >> > If you do not work on an important problem, it's unlikely >> > you'll do important work. -- R. Hamming >> > http://www.cs.virginia.edu/~robins/YouAndYourResearch.html >> > >> > -- >> > You received this message because you are subscribed to the Google >> Groups "sympy" group. >> > To unsubscribe from this group and stop receiving emails from it, send >> an email to sympy+unsubscr...@googlegroups.com. >> > To post to this group, send email to sympy@googlegroups.com. >> > Visit this group at http://groups.google.com/group/sympy. >> > To view this discussion on the web visit >> https://groups.google.com/d/msgid/sympy/20150216095625.GQ6082%40machine.or.cz >> . >> > For more options, visit https://groups.google.com/d/optout. >> >> -- >> You received this message because you are subscribed to the Google Groups >> "sympy" group. >> To unsubscribe from this group and stop receiving emails from it, send an >> email to sympy+unsubscr...@googlegroups.com. >> To post to this group, send email to sympy@googlegroups.com. >> Visit this group at http://groups.google.com/group/sympy. >> To view this discussion on the web visit >> https://groups.google.com/d/msgid/sympy/CADDwiVDaF6Q%2BO1RTJRRiEfDiG23bJicPZFCXwq3CQXJfys2XZw%40mail.gmail.com >> . >> For more options, visit https://groups.google.com/d/optout. >> > > -- You received this message because you are subscribed to the Google Groups "sympy" group. To unsubscribe from this group and stop receiving emails from it, send an email to sympy+unsubscr...@googlegroups.com. To post to this group, send email to sympy@googlegroups.com. Visit this group at http://groups.google.com/group/sympy. To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/CAP7f1AgGfqnqxAMinfMZquVqxKy4j8i_7CXSeVJ%2BxWYuouGK-A%40mail.gmail.com. For more options, visit https://groups.google.com/d/optout.
