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/
-~----------~----~----~----~------~----~------~--~---

Reply via email to