Raymond Toy | 1 Sep 14:46 2006
Picon

Re: adding floats. a peculiarity, and a proposal

>>>>> "Richard" == Richard Fateman <fateman <at> cs.berkeley.edu> writes:

    Richard> in 5.9.3, xmaxima
    Richard> (%i8) sum(0.1d0,i,1,1000);
    Richard> (%o8) 			       99.99999999999859
    Richard> (%i9) sum(0.1d0,i,1,1000000);
    Richard> (%o9) 			       100000.0000013329

    Richard> It's like that old joke "we lose a little on each sale, but we make it up on
    Richard> the volume..."

Haha!

    Richard> I have (yet another) proposal for floats for people to think
    Richard> about, courtesy of W. Kahan..

    Richard> In one sentence:

    Richard>  ALL user floating-point format input be stored as infinitely precise
    Richard> rational numbers until
    Richard> there is a first conversion to something that is used computationally, at
    Richard> which point it is
    Richard> converted to the appropriate type, precision, etc.

So in the examples above, what would be the result? 100 and 1000000?
100d0 and 100000d0?

Not sure I like this.  If I said 0.1d0, then I want 0.1d0.  If I
wanted an exact rational I would have said 1/10.  However, it would be
useful to have, say, 0.1234567r0 mean the exact rational
(Continue reading)

Stavros Macrakis | 1 Sep 16:19 2006
Picon

Re: adding floats. a peculiarity, and a proposal

>  ALL user floating-point format input be stored as infinitely precise
> rational numbers until there is a first conversion to something
> that is used computationally, at which point it is converted to the
> appropriate type, precision, etc.

I am not sure exactly what you are proposing.  If 0.1b0 is stored as
1/10 but converted when it is "used computationally", then you will
get the same results as above, but with a large overhead:
sum(0.1b0,i,1,10) would have to convert 1/10 to binary form 10 times.

And what is the "appropriate type"?  The value of fpprec at the time
of conversion?  That is the current semantics for converting exact
numbers to bfloats.

Perhaps the proposal is that bfloat constants be represented as
decimal floating-point?  That certainly has advantages (and
disadvantages).  And I believe that bfloat already supports base-10
internally in some cases at least.

As for the specific examples of sum, I wonder if they would come out
differently if we fixed the known bugs in bfloat?  There are
documented cases where bfloat does not round correctly.

By the way, I was curious to see what would happen to computation
times if I used exact (rational) arithmetic (sum(1/i,i,1,N)) vs.
bfloat (sum(1.0b0/i,i,1,N)).  I expected that with the growth in
denominator size, rational would quickly become much more expensive
than bfloat (default fpprec=16).  In fact, the crossover is at about
N=10000.

(Continue reading)

Richard Fateman | 1 Sep 16:36 2006
Picon

Re: adding floats. a peculiarity, and a proposal


> -----Original Message-----
> From: Stavros Macrakis [mailto:macrakis <at> gmail.com] 
> Sent: Friday, September 01, 2006 7:20 AM
> To: Richard Fateman
> Cc: maxima
> Subject: Re: [Maxima] adding floats. a peculiarity, and a proposal
> 
> >  ALL user floating-point format input be stored as 
> infinitely precise 
> > rational numbers until there is a first conversion to 
> something that 
> > is used computationally, at which point it is converted to the 
> > appropriate type, precision, etc.
> 
> I am not sure exactly what you are proposing.  If 0.1b0 is 
> stored as 1/10 but converted when it is "used 
> computationally", then you will get the same results as 
> above, but with a large overhead:
> sum(0.1b0,i,1,10) would have to convert 1/10 to binary form 10 times.

There are areas that would have to be figured out, maybe heuristically

0.1b0 with default fpprec is ((BIGFLOAT SIMP 56) 57646075230342349 -3) 
My proposal would change that to either ((BIGFLOAT SIMP 56 1/10)
57646075230342349 -3) 
or maybe ((rat simp ((BIGFLOAT SIMP 56) 57646075230342349 -3) )  1 10)

if you added two of these together, it could be done as exact rationals or
as bigfloats.
(Continue reading)

Richard Fateman | 1 Sep 16:44 2006
Picon

Re: adding floats. a peculiarity, and a proposal


> -----Original Message-----
> From: Raymond Toy [mailto:raymond.toy <at> ericsson.com] 
....> 
> Not sure I like this.  If I said 0.1d0, then I want 0.1d0.  
> If I wanted an exact rational I would have said 1/10.  

If you have an explicit marker "I want a double" that would make sense.

If you have a somewhat vaguer marker "I want to use decimal positional
notation,
(perhaps with an exponent e)" then what does that mean?
Does 0.1  or 0.1e0 mean  1/10 or does it mean
0.10000000000000000034694469519536141888238489627838134765625..... ?

I think that 1/10 is less confusing for the computationally naive; the
expert who
knows about floating point arithmetic can be told to use d.  The person who
knows
only a little about floats is often confused anyway :)

> However, it would be useful to have, say, 0.1234567r0 mean 
> the exact rational 1234567/100000000, without me having to 
> write the denominator explicitly.

This is possible too:  how many exponent markers do we want? eEdDsSbBrR?

RJF
Alasdair McAndrew | 2 Sep 00:25 2006
Picon

Mod, ?mod and matrices

I have a little annoyance with the mod function.  As you know mod(x,p):

      is defined as an integer <k> in
     `[-(<p>-1)/2, ..., 0, ..., (<p>-1)/2]' when <p> is odd, or
     `[-(<p>/2 - 1), ..., 0, ...., <p>/2]' when <p> is even, such that
     `<a> <p> + <k>' equals <n> for some integer <a>

This seems to me to be non-standard; more generally mod(x,n) is defined do be in the range 0.. n-1.  Now, I can use ?mod to do this, except that ?mod doesn't apply to matrices in the same way mod does:

--
(%i19) M:genmatrix(lambda([i,j],random(100)),3,3);
                                [ 60  71  12 ]
                                [            ]
(%o19)                          [ 9   29  18 ]
                                [            ]
                                [ 26  95  53 ]
(%i20) mod(M,10);
                               [  0    1    2  ]
                               [               ]
(%o20)                         [ - 1  - 1  - 2 ]
                               [               ]
                               [ - 4   5    3  ]
(%i21) ?mod(M, 10);
Maxima encountered a Lisp error:


MOD: (($MATRIX SIMP) ((MLIST SIMP) 60 71 12) ((MLIST SIMP) 9 29 18)
 ((MLIST SIMP) 26 95 53)) is not a REAL

Automatically continuing.
To reenable the Lisp debugger set *debugger-hook* to nil.
(%i22)  matrixmap(lambda([x],?mod(x,10)),M);
                                  [ 0  1  2 ]
                                  [         ]
(%o22)                            [ 9  9  8 ]
                                  [         ]
                                  [ 6  5  3 ]
--
Is there any way of coercing map or ?map to produce the last answer without having to fiddle with matrixmap and lambda definitions?

Thanks,
Alasdair







_______________________________________________
Maxima mailing list
Maxima <at> math.utexas.edu
http://www.math.utexas.edu/mailman/listinfo/maxima
Alasdair McAndrew | 2 Sep 02:33 2006
Picon

Re: Mod, ?mod and matrices

Sorry: in the last line I meant mod or ?mod, not map or ?map.

-A.

_______________________________________________
Maxima mailing list
Maxima <at> math.utexas.edu
http://www.math.utexas.edu/mailman/listinfo/maxima
Barton Willis | 2 Sep 02:59 2006

Re: Mod, ?mod and matrices

-----maxima-bounces <at> math.utexas.edu wrote: -----

>To: "maxima list" <maxima <at> math.utexas.edu>
>From: "Alasdair McAndrew" <amca01 <at> gmail.com>

>I have a little annoyance with the mod
>function.  As you know mod(x,p):
>
>      is defined as an integer <k> in
>
>     `[-(-1)/2, ..., 0, ..., (-1)/2]' when
><p> is odd, or
>
>     `[-(/2 - 1), ..., 0, ...., /2]' when
><p> is even, such that

This is no longer the case (starting with maybe version 5.9.3 or so):

 -- Function: mod (<x>, <y>)
     If <x> and <y> are real numbers and <y> is nonzero, return `<x> -
     <y> * floor(<x> / <y>)'.  Further for all real <x>, we have `mod
     (<x>, 0) = <x>'. For a discussion of the definition `mod (<x>, 0)
     = <x>', see Section 3.4, of "Concrete Mathematics," by Graham,
     Knuth, and Patashnik.

To apply mod to each element of a matrix, you'll still need to use
matrixmap
and a lambda form. Doing mod(matrix([...]),10) makes a mess...

Barton
Picon

Switching from ode2.lisp to ode2.mac

The ode2() function returns Bessel function solutions of form 
%j[nu](x), but the function %j has recently been moved from 
maxima.

Rather than modify the translated file ode2.lisp, I'd prefer to
switch to the original maxima file ode2.mac and modify that.
I've tested this change on the maxima and share/contrib/diffequation
testsuites.  The results are identical and the the timings are
similar.

Thoughts?

NOTICE
This e-mail and any attachments are private and confidential and may contain privileged information. If
you are not an authorised recipient, the copying or distribution of this e-mail and any attachments is
prohibited and you must not read, print or act in reliance on this e-mail or attachments.
This notice should not be removed.
Stavros Macrakis | 2 Sep 15:05 2006
Picon

Re: Mod, ?mod and matrices

On 9/1/06, Barton Willis <willisb <at> unk.edu> wrote:
> To apply mod to each element of a matrix, you'll still need to use
> matrixmap and a lambda form. Doing mod(matrix([...]),10) makes a mess...

Well, that's a bug we should fix... both mod and floor should
distribute over matrix and list.  The current result is silly.
Robert Dodier | 2 Sep 17:19 2006
Picon

Re: Switching from ode2.lisp to ode2.mac

On 9/2/06, Billinghurst, David (CALCRTS)
<david.billinghurst <at> comalco.riotinto.com.au> wrote:

> Rather than modify the translated file ode2.lisp, I'd prefer to
> switch to the original maxima file ode2.mac and modify that.
> I've tested this change on the maxima and share/contrib/diffequation
> testsuites.  The results are identical and the the timings are
> similar.

Makes sense to me. Seems like the user-level functions
in ode2.mac should be given an autoload property --
otherwise it would be necessary to load the script
explicitly to get the functions.

best,
Robert

Gmane