On Dec 3, 2007 8:42 AM, Jason Grout <[EMAIL PROTECTED]> wrot > > > William Stein wrote: > > On Dec 3, 2007 8:13 AM, Bill Hart <[EMAIL PROTECTED]> wrote: > >> I did try to check that Mathematica was getting the right answer, but > >> I had no luck. I don't know how to convert a mathematica matrix into > >> ordinary matrix form in SAGE, so when I do the comparison it always > >> just says false. > > > > Damn it -- you're right. This isn't the first time I've been bitten by > > this. Mathematica is just doing the componentwise product! > > > > sage: mm = mathematica(n) > > sage: n > > [ 1/20 -1/10 1/20] > > [ -1/10 -37/15 47/30] > > [ 1/20 67/30 -77/60] > > sage: mm^2 > > {{1/400}, {1/100}, {1/400}, {1/100}, {1369/225}, {2209/900}, {1/400}, > > {4489/900}, {5929/3600}} > > sage: n^2 > > [ 3/200 53/150 -131/600] > > [ 8/25 1439/150 -147/25] > > [ -57/200 -419/50 3089/600] > > > > > > Why the frick do so many mathematical software systems define > > A^n for A a matrix to be componentwise powering. Other programs, > > e.g., numpy define A*B on numpy arrays to be componentwise > > multiplication. > > > > > > So -- does anybody know how to raise a matrix to a power in > > Mathematica, since I sure don't. > > > > Mathematica treats matrices as lists of lists. Since the product of a > list and another list is the component-wise product, the "*" means > hadamard multiplication. > > You probably want the "." operator, which does matrix multiplication: > > > {{1, 0}, {0, 1}}.{{1, 2}, {3, 4}} > > > {{1, 2}, {3, 4}} > > or Dot[{{1,0},{0,1}}, {{1,2},{3,4}}] > > Yes, this has bitten me several times too. You could also use: > > MatrixPower[{{1,2,3},{4,5,6},{7,8,9}}, 20000] > > to compute the matrix product power. For benchmarking, this is probably > the function to use.
Thanks. I first found that there was a bug in the sage --> mathematica matrix conversion code: http://trac.sagemath.org/sage_trac/ticket/1382 Then I tried using what you said above directly to get timings: sage: m = mathematica('{{1/20, -1/10, 1/20}, {-1/10, -37/15, 47/30}, {1/20, 67/30,-77/60}}') sage: m.name() 'sage5' sage: time a = mathematica('MatrixPower[sage5, 20000]') CPU times: user 0.00 s, sys: 0.00 s, total: 0.00 s Wall time: 0.81 sage: time a = mathematica('MatrixPower[sage5, 200000]') CPU times: user 0.00 s, sys: 0.00 s, total: 0.00 s Wall time: 38.94 Recall that in Sage (n is the same matrix as sage5 above): sage: time a = n^20000 CPU times: user 0.16 s, sys: 0.00 s, total: 0.16 s Wall time: 0.16 sage: time a = n^200000 CPU times: user 4.38 s, sys: 0.05 s, total: 4.43 s Wall time: 4.42 So Sage is *already* the fastest in the world at this benchmark by far. (OK, so I didn't test maple, since I don't know how, but I doubt it is better...) -- William --~--~---------~--~----~------------~-------~--~----~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to [EMAIL PROTECTED] For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://sage.scipy.org/sage/ and http://modular.math.washington.edu/sage/ -~----------~----~----~----~------~----~------~--~---