Just to be more clear, I'm not asking about powering float Matrices, per se.

Other plausible techniques are known. Eg. since "almost all" float Matrices are diagonalizable, the product (P^(-1).D.P)^n telescopes down to P^(-1).D^n.P and  only  the scalar  eigenvalues  would be powered.

I'm wondering more about the prospects for optimization following enhancement for iterative specialization and, perhaps more difficult, handling the in-place semantic.

acer

I thoroughly enjoyed reading the thesis, and I look forward to poking at the code.

I have just one or two quite specific questions, which I'd like to ask now if I may. In subsection 7.2 there is an example on binary powering (Listings 7.7 and 7.8). In the specialized code of Listing 7.8 there are more local variables than are strictly necessary. Ok, so locals are not too expensive, in themselves. But what about optimizing for the least number of such locals? Is that something that might fall under the possible future work of "Iterative Specialization" as mentioned in subsection 8.3.2? I ask because I'm thinking of the case of powering a Matrix rather than an scalar expression. In the hardware datatype case, the entries of Matrices are never garbage, but their container Matrices are. So the more locals that get produced (and assigned to, rather than copied to in-place), the more expensive it is.

I'll try to explain by example. Consider the binary powering of hardware (any, actually) Matrices as done in Maple 9.5. It is akin to Listing 7.7, and it produces whole new Matrices galore. And it's really slow, on account of all the ensuing garbage collection.

    |\^/|     Maple 9.5 (IBM INTEL LINUX)
._|\|   |/|_. Copyright (c) Maplesoft, a division of Waterloo Maple Inc. 2004
 \  MAPLE  /  All rights reserved. Maple is a trademark of
 <____ ____>  Waterloo Maple Inc.
      |       Type ? for help.
> interface(verboseproc=3):
> eval(`rtable/Power/Power`);
proc(M, n)
local P;
option `Copyright (c) 1999 Waterloo Maple Inc. All rights reserved.`;
    if n = 1 then M
    elif n = 2 then LinearAlgebra:-MatrixMatrixMultiply(M, M)
    else
        P := `rtable/Power/Power`(M, iquo(n, 2));
        P := LinearAlgebra:-MatrixMatrixMultiply(P, P);
        if n mod 2 = 1 then
            LinearAlgebra:-MatrixMatrixMultiply(P, M, 'inplace' = true)
        end if;
        P
    end if
end proc

Now look at the current implementation. I won't show the whole listing, but it can be described easily enough.

interface(verboseproc=3);
eval(`rtable/Power/float`);

It has distinct "specializations" for powers 1,2,3,4, and 5. And each of those minimizes the number of allocated hardware Matrices. Matrices are never assigned to locals except when first created. No other instances of assigning a Matrix to a variable occurs, and all updates are done by in-place copying. And then it has a version for powers 6 and higher, which follows a pattern (common for any p>6) and thus uses a common scheme minimized for the total number of allocated Matrices.

What I'm wondering is, is there any hope of generating programs such as that using MapleMix techniques (present, or with any future idea from section 8)? I would like to create such optimized programs, but I don't want to have to work it out by hand.

acer

Just a minor point: that exact syntax has been available for many releases now.

What's new in Maple 12 (I think) is this,

> A := Array(1..10,1..10,1..10, (i,j,k)->i+j/k ):

> A[1,1,1..];
               [                                            11]
               [2, 3/2, 4/3, 5/4, 6/5, 7/6, 8/7, 9/8, 10/9, --]
               [                                            10]

> A[1,1,4..]:

> Vector[row](A[1,1,4..]);
                      [                               11]
                      [5/4, 6/5, 7/6, 8/7, 9/8, 10/9, --]
                      [                               10]

acer

Just a minor point: that exact syntax has been available for many releases now.

What's new in Maple 12 (I think) is this,

> A := Array(1..10,1..10,1..10, (i,j,k)->i+j/k ):

> A[1,1,1..];
               [                                            11]
               [2, 3/2, 4/3, 5/4, 6/5, 7/6, 8/7, 9/8, 10/9, --]
               [                                            10]

> A[1,1,4..]:

> Vector[row](A[1,1,4..]);
                      [                               11]
                      [5/4, 6/5, 7/6, 8/7, 9/8, 10/9, --]
                      [                               10]

acer

The try..catch mechanism allows one to control the program flow, even in the presence of an error.

For example,

> p := proc(x)
> local z;
> try
>   if 0 < x then error "whoops";
>   else z := x^2;
>   end if;
> catch "whoops":
>   print("found a whoops");
>   z := sin(z);
> end try;
> end proc:

> p(-4);
                                      16
 
> p(4);
                               "found a whoops"
 
                                    sin(z)

So, simply catch the particular error message which LPSolve issues when it's infeasible.

See the ?try help-page.

acer

The try..catch mechanism allows one to control the program flow, even in the presence of an error.

For example,

> p := proc(x)
> local z;
> try
>   if 0 < x then error "whoops";
>   else z := x^2;
>   end if;
> catch "whoops":
>   print("found a whoops");
>   z := sin(z);
> end try;
> end proc:

> p(-4);
                                      16
 
> p(4);
                               "found a whoops"
 
                                    sin(z)

So, simply catch the particular error message which LPSolve issues when it's infeasible.

See the ?try help-page.

acer

Am I just missing it, or is there no reference or link to ?Rounding from either the ?evalf or ?evalf,details help-pages? Also, maybe the ?Rounding page could benefit from an example such as,

> evalf[5](Pi);
                                    3.1416

> Rounding:=0:
> forget(evalf):
> evalf[5](Pi);
                                    3.1415

Yes, as Alec points out iquo (or irem) is usual. See also,

interface(verboseproc=3):
eval(`convert/base`);

acer

Am I just missing it, or is there no reference or link to ?Rounding from either the ?evalf or ?evalf,details help-pages? Also, maybe the ?Rounding page could benefit from an example such as,

> evalf[5](Pi);
                                    3.1416

> Rounding:=0:
> forget(evalf):
> evalf[5](Pi);
                                    3.1415

Yes, as Alec points out iquo (or irem) is usual. See also,

interface(verboseproc=3):
eval(`convert/base`);

acer

That's a good point, about the timing function. A more convenient version of Maple's time() command, which returns both the result and the elapsed time, would be useful.

As far as the speed goes, does anyone know how many pre-computed digits of Pi Mathematica has stored (hard-coded)? It's even more interesting if one knows that both systems had to compute the answer in full.

acer

That's a good point, about the timing function. A more convenient version of Maple's time() command, which returns both the result and the elapsed time, would be useful.

As far as the speed goes, does anyone know how many pre-computed digits of Pi Mathematica has stored (hard-coded)? It's even more interesting if one knows that both systems had to compute the answer in full.

acer

Nice.

My rather facile answer above was due to my mistaking it as a Student forum question, possibly an exercise from calculus. Sorry about that.

acer

Nice.

My rather facile answer above was due to my mistaking it as a Student forum question, possibly an exercise from calculus. Sorry about that.

acer

Well, it's easy to be fast when looking up precomputed results.

Fast lookup, even if it requires a base conversion, is about the opposite end of the interest spectrum for me as is the BBP research (or similar, where most earlier digits aren't even computed when searching for some high nth).

acer

Well, it's easy to be fast when looking up precomputed results.

Fast lookup, even if it requires a base conversion, is about the opposite end of the interest spectrum for me as is the BBP research (or similar, where most earlier digits aren't even computed when searching for some high nth).

acer

The default for Standard is a thing called a Document. But there is also a Worksheet.

One can always open a new one of either, using the File->Open menubar item.

I set my global preference to Worksheet (again, somewhere under Tools->Options).

One can use either 2D or 1D input in a Worksheet.

So, really, there are two sorts of option/choice for the Standard GUI. Both can be set separately. Bot can be set as a preference (globally, for all new sessions).

• 1D Maple input vs 2D Math input
• Worksheet vs Document

Writing cool applications with embedded components linked to live(!) tyepset 2D math equations is where I might use a Document and 2D Math input. But for everyday use, and especially for coding, I use a Worksheet and 1D Maple input.

I like 1D Maple input because what I see is what is actually there, unambiguously.

acer