On Jun 6, 2010, at 6:29 AM, Mateusz Paprocki wrote:

> Hi,
> 
> On Sat, Jun 05, 2010 at 06:10:51PM -0600, Aaron S. Meurer wrote:
>> I have rewritten heurisch() using Poly, but it is now considerably slower.  
>> Here is the reason (pardon the long expression).  This is take from 
>> debugging integrate(exp(x)*sin(x)*x**2) in my branch:
>> 
>> In [57]: poly_part
>> Out[57]: Poly(_x0**5*_A79 + _x0**4*_x1*_A53 + _x0**4*_x2*_A105 + 
>> _x0**4*_x3*_A40 + _x0**4*_A57 + _x0**3*_x1**2*_A58 + _x0**3*_x1*_x2*_A116 + 
>> _x0**3*_x1*_x3*_A106 + _x0**3*_x1*_A77 + _x0**3*_x2**2*_A14 + 
>> _x0**3*_x2*_x3*_A96 + _x0**3*_x2*_A75 + _x0**3*_x3**2*_A121 + 
>> _x0**3*_x3*_A50 + _x0**3*_A59 + _x0**2*_x1**3*_A13 + _x0**2*_x1**2*_x2*_A110 
>> + _x0**2*_x1**2*_x3*_A82 + _x0**2*_x1**2*_A30 + _x0**2*_x1*_x2**2*_A93 + 
>> _x0**2*_x1*_x2*_x3*_A112 + _x0**2*_x1*_x2*_A109 + _x0**2*_x1*_x3**2*_A89 + 
>> _x0**2*_x1*_x3*_A123 + _x0**2*_x1*_A10 + _x0**2*_x2**3*_A4 + 
>> _x0**2*_x2**2*_x3*_A15 + _x0**2*_x2**2*_A7 + _x0**2*_x2*_x3**2*_A27 + 
>> _x0**2*_x2*_x3*_A18 + _x0**2*_x2*_A72 + _x0**2*_x3**3*_A68 + 
>> _x0**2*_x3**2*_A70 + _x0**2*_x3*_A0 + _x0**2*_A11 + _x0*_x1**4*_A73 + 
>> _x0*_x1**3*_x2*_A32 + _x0*_x1**3*_x3*_A43 + _x0*_x1**3*_A88 + 
>> _x0*_x1**2*_x2**2*_A26 + _x0*_x1**2*_x2*_x3*_A104 + _x0*_x1**2*_x2*_A23 + 
>> _x0*_x1**2*_x3**2*_A60 + _x0*_x1**2*_x3*_A35 + _x0*_x1**2*_A8 + 
>> _x0*_x1*_x2**3*_A103 + _x0*_x1*_x2**2*_x3*_A122 + _x0*_x1*_x2**2*_A76 +
>> _x0*_x1*_x2*_x3**2*_A117 + _x0*_x1*_x2*_x3*_A54 + _x0*_x1*_x2*_A12 + 
>> _x0*_x1*_x3**3*_A69 + _x0*_x1*_x3**2*_A114 + _x0*_x1*_x3*_A37 + _x0*_x1*_A74 
>> + _x0*_x2**4*_A95 + _x0*_x2**3*_x3*_A80 + _x0*_x2**3*_A16 + 
>> _x0*_x2**2*_x3**2*_A38 + _x0*_x2**2*_x3*_A39 + _x0*_x2**2*_A3 + 
>> _x0*_x2*_x3**3*_A101 + _x0*_x2*_x3**2*_A67 + _x0*_x2*_x3*_A85 + _x0*_x2*_A56 
>> + _x0*_x3**4*_A45 + _x0*_x3**3*_A21 + _x0*_x3**2*_A98 + _x0*_x3*_A108 + 
>> _x0*_A125 + _x1**5*_A102 + _x1**4*_x2*_A51 + _x1**4*_x3*_A6 + _x1**4*_A25 + 
>> _x1**3*_x2**2*_A36 + _x1**3*_x2*_x3*_A124 + _x1**3*_x2*_A47 + 
>> _x1**3*_x3**2*_A33 + _x1**3*_x3*_A84 + _x1**3*_A81 + _x1**2*_x2**3*_A91 + 
>> _x1**2*_x2**2*_x3*_A94 + _x1**2*_x2**2*_A29 + _x1**2*_x2*_x3**2*_A107 + 
>> _x1**2*_x2*_x3*_A92 + _x1**2*_x2*_A20 + _x1**2*_x3**3*_A87 + 
>> _x1**2*_x3**2*_A31 + _x1**2*_x3*_A52 + _x1**2*_A115 + _x1*_x2**4*_A24 + 
>> _x1*_x2**3*_x3*_A22 + _x1*_x2**3*_A55 + _x1*_x2**2*_x3**2*_A64 + 
>> _x1*_x2**2*_x3*_A100 + _x1*_x2**2*_A48 + _x1*_x2*_x3**3*_A119 + 
>> _x1*_x2*_x3**2*_A42 + _x1*_x2*_x3*_A28 + _x1*_x2*_A49 + _
>> x1*_x3**4*_A1 + _x1*_x3**3*_A62 + _x1*_x3**2*_A41 + _x1*_x3*_A113 + _x1*_A34 
>> + _x2**5*_A44 + _x2**4*_x3*_A99 + _x2**4*_A83 + _x2**3*_x3**2*_A63 + 
>> _x2**3*_x3*_A111 + _x2**3*_A78 + _x2**2*_x3**3*_A5 + _x2**2*_x3**2*_A2 + 
>> _x2**2*_x3*_A97 + _x2**2*_A120 + _x2*_x3**4*_A65 + _x2*_x3**3*_A66 + 
>> _x2*_x3**2*_A46 + _x2*_x3*_A118 + _x2*_A71 + _x3**5*_A19 + _x3**4*_A90 + 
>> _x3**3*_A61 + _x3**2*_A17 + _x3*_A86 + _A9, _x0, _x1, _x2, _x3, _A0, _A1, 
>> _A2, _A3, _A4, _A5, _A6, _A7, _A8, _A9, _A10, _A11, _A12, _A13, _A14, _A15, 
>> _A16, _A17, _A18, _A19, _A20, _A21, _A22, _A23, _A24, _A25, _A26, _A27, 
>> _A28, _A29, _A30, _A31, _A32, _A33, _A34, _A35, _A36, _A37, _A38, _A39, 
>> _A40, _A41, _A42, _A43, _A44, _A45, _A46, _A47, _A48, _A49, _A50, _A51, 
>> _A52, _A53, _A54, _A55, _A56, _A57, _A58, _A59, _A60, _A61, _A62, _A63, 
>> _A64, _A65, _A66, _A67, _A68, _A69, _A70, _A71, _A72, _A73, _A74, _A75, 
>> _A76, _A77, _A78, _A79, _A80, _A81, _A82, _A83, _A84, _A85, _A86, _A87, 
>> _A88, _A89, _A90, _A91, _A92, _A93, _A94, _A95, _A96, _A97, _A98, _A99, 
>> _A100, _
>> A101, _A102, _A103, _A104, _A105, _A106, _A107, _A108, _A109, _A110, _A111, 
>> _A112, _A113, _A114, _A115, _A116, _A117, _A118, _A119, _A120, _A121, _A122, 
>> _A123, _A124, _A125, domain='ZZ')
>> 
>> In [58]: b
>> Out[58]: Poly(_x0, _x0, _x1, _x2, _x3, _A0, _A1, _A2, _A3, _A4, _A5, _A6, 
>> _A7, _A8, _A9, _A10, _A11, _A12, _A13, _A14, _A15, _A16, _A17, _A18, _A19, 
>> _A20, _A21, _A22, _A23, _A24, _A25, _A26, _A27, _A28, _A29, _A30, _A31, 
>> _A32, _A33, _A34, _A35, _A36, _A37, _A38, _A39, _A40, _A41, _A42, _A43, 
>> _A44, _A45, _A46, _A47, _A48, _A49, _A50, _A51, _A52, _A53, _A54, _A55, 
>> _A56, _A57, _A58, _A59, _A60, _A61, _A62, _A63, _A64, _A65, _A66, _A67, 
>> _A68, _A69, _A70, _A71, _A72, _A73, _A74, _A75, _A76, _A77, _A78, _A79, 
>> _A80, _A81, _A82, _A83, _A84, _A85, _A86, _A87, _A88, _A89, _A90, _A91, 
>> _A92, _A93, _A94, _A95, _A96, _A97, _A98, _A99, _A100, _A101, _A102, _A103, 
>> _A104, _A105, _A106, _A107, _A108, _A109, _A110, _A111, _A112, _A113, _A114, 
>> _A115, _A116, _A117, _A118, _A119, _A120, _A121, _A122, _A123, _A124, _A125, 
>> domain='ZZ[_x0]')
>> 
>> In [59]: %timeit poly_part*b
>> 1 loops, best of 3: 2.02 s per loop
>> 
>> In [60]: %timeit Poly((poly_part.as_basic()*b.as_basic()).expand())
>> 1 loops, best of 3: 870 ms per loop
>> 
>> So it is much faster to convert the Polys to basic, expand them using 
>> inefficient expand(), then put them back into a poly.  There has to be a 
>> better way to do this.  I know that in the long term, sdp_mul will be the 
>> solution (for monomial multiplication like this, it is almost 
>> instantaneous).  Is there a way to make dmp multiplication faster?  And if 
>> not, what exactly are the plans for sdp, i.e., assuming we had a working SDP 
>> class, would it completely replace DMP, or was the plan for Poly to use 
>> which ever one was more efficient for whichever particular operation?  And 
>> regarding my GSoC project, would it be better for me to skip this for now 
>> and start work on the full algorithm, or for me to do some work on sdp to 
>> make things faster?
>> 
> 
> My suggestion is to focus more on the algorithm(s) and less on the
> technical details at this point. We know that DMP representation is
> a bad choice for polynomials with so many variables (there are several
> issues open concerning this). My plan is to employ SDP as the default
> representation and leave DMP as backup (there are cases in which DMP
> can be several times faster than SDP). Basically, when are implementing
> so high-level algorithms, as Risch algorithm is, you shouldn't be very
> concerned what is happening on the lower levels. I'm not sure when
> exactly I will have enough time to finalize SDP code, but it should
> happen before August. If you will find time to work on this or other
> issues during GSoC, that's fine, but I would like to see the algorithms
> working --- speed optimisations can be done later.

OK.  That is what I figured the answer would be.  It's best to focus on 
correctness first, then speed, digging into the lower levels only if I need to 
implement a feature or fix a bug.  Don't worry, though; time hasn't been 
wasted, because now I am quite familiar with the flow of the heurisch code, 
much of which will be very similar if not identical to code I will need for the 
full algorithm.    

Actually, normally I wouldn't bother too much with the speed at this point, but 
the Polys are so unbearably slow in this case.  So I guess I will just keep the 
Polified heurisch on the back burner somewhere and see if it performs better 
with SDP Polys in the future.  

Aaron Meurer
> 
>> Sorry, these are mostly directed toward Mateusz, but if anyone else wants to 
>> chime in, feel free.
>> 
>> Aaron Meurer   
>> 
>> -- 
>> 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.
>> 
> 
> -- 
> Mateusz
> 

-- 
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.

Reply via email to