🚀 go-pugleaf

Glen Herrmannsfeldt wrote:

> As far as I know, Fortran still requires floating point multiply to
> generate a result with no more precision than the arguments that went into
> it. W.Kahan describes well the advantage of summing a double precision
> product of single precision arguments, yet the Fortran standard doesn't
> seem to
> follow it.  (Someone with a more recent standard can comment on this.)
> 
I doubt that all the machines we've been using violate the Fortran 
standard.  Didn't f77 introduce DPROD() specifically to require a double 
precision product, in case that should be required?
> If someone wanted to comment about other hardware or compilers in use in
> the 1960's and 1970's, that might also be useful.
> 
>
The GE600/Honeywell6000, and probably the Univac 36-bit machines, generated 
an 80-bit product.

-- 
Tim Prince

In article <dG4db.598668$uu5.97984@sccrnsc04>,
"Glen Herrmannsfeldt" <gah@ugcs.caltech.edu> wrote:

 > An article by W. Kahan titled "How Java's Floating-Point Hurts Everyone
 > Everywhere" (and which describes much more than Java)  includes a
discussion
 > about the results of floating point multiply.

Where might one find this article?  Online?

 > According to his description, until 1963 IBM Fortran compilers would
 > generate double precision results for floating point multiply of two
single
 > precision numbers, at least when it was being summed into a double
precision
 > variable.  He then explains that IBM changed this in IBSYS 13,
anticipating
 > S/360.

They didn't change the way the 7090 and 7094 computed the result, which was
that the LSB of the 36x36 multiply was physically in a register --
available.  If the programmer had specified that the result of
single*single was to be placed in a double, then it could be, with no
further work.  Obviously, single=single*single didn't produce a double
result.

 > As far as I can tell, the Fortran-66 standard requires a single
precision
 > product when multiplying single precision numbers, though I don't know
how
 > much impact IBM had on Fortran-66.

I'd appreciate an exact quote on that.  I used to know X3.9-1966 awfully
well, and I don't recall that the Honeywell gear that I used after the 7094
and 360 and CDC 3300 had a standards violation in that regard.  My copy is
filed away, though.

 > The reason why this comment is surprising to me is that S/360 (along
with
 > S/370 and ESA/390) do not have a single precision multipy instruction
with a
 > single precision product.  The product is always double precision.  I
have
 > checked the code generated by both the G and H compilers, and both,
after
 > generating the multiply instruction generate at least two more
instructions
 > to truncate the product to single precision.

I wonder why they bother.  If the result is to be stored, the extra
precision goes away.  If it's to continue in a single-precision
computation, it doesn't help much -- only the first bit of the LSBs has any
effect on the next multiplication, and rounding introduces a possible error
of the same magnitude.

 > As far as I know, Fortran still requires floating point multiply to
generate
 > a result with no more precision than the arguments that went into it.
 > W.Kahan describes well the advantage of summing a double precision
product
 > of single precision arguments, yet the Fortran standard doesn't seem to
 > follow it.  (Someone with a more recent standard can comment on this.)

Honeywell had what are called overlength registers, as Tim Prince has
mentioned.  A 36-bit float was actually an 8-bit exponent and a 28-bit
fraction.  On loading these into registers, they went into an 8-bit
E-register, and a 36-bit (machine word-length) register.  A following
multiplication, for example, produced a normalized result by updating the
Exponent register, with overlength arithmetic product in the full-length
register.  That overlength product or other intermediate result was
maintained so long as the computation did not send a result to storage, but
a floating-store rounded the value properly and put it back into a 36-bit
(8+28 bit) float.  It was a terrific system and produced very good, very
acceptable to programmers, results.  I'm frankly skeptical that this
violated X3.9-1966.

 > I would, in the end, like to thank W.Kahan for his unceasing drive to
 > improve the results of floating point hardware and software, and hope
for
 > only constructive discussion of this problem.

That cite, so the rest of us could look at the article to which you're
referring, would be real nice.

I had to respond to government (US national government) queries about the
Honeywell overlength registers on a regular basis, and I used to include
with my replies to specific questions a copy of a paper by W. D. Maurer
(I'm pretty sure that's the one), from Communications of the ACM.  His
paper was on the results of floating point arithmetic operations on many
platforms, with the overlength issue a minor mention near the end of the
paper.  I don't have a precise cite but the paper came out in roughly 1972.
 It may just have been my magic personality or my renowned writing skills,
but I like to think that we pretty well settled all the concerns they had
about Honeywell's somewhat unusual arithmetic arrangement.  We were never
disqualified from a bid, and I don't recall any case in which we were
either marked-down or questioned further, once the potential customer had a
chance to review our specific comments and the paper from CACM.

      Edward Rice

"Edward Rice" <ehrice@his.com> wrote in message
news:BB9A8C689668106A5D@10.0.1.2...
> In article <dG4db.598668$uu5.97984@sccrnsc04>,
> "Glen Herrmannsfeldt" <gah@ugcs.caltech.edu> wrote:
>
>  > An article by W. Kahan titled "How Java's Floating-Point Hurts Everyone
>  > Everywhere" (and which describes much more than Java)  includes a
> discussion
>  > about the results of floating point multiply.
>
> Where might one find this article?  Online?

A google search with the title should find it, but...

http://www.cs.berkeley.edu/~wkahan/JAVAhurt.pdf

>  > According to his description, until 1963 IBM Fortran compilers would
>  > generate double precision results for floating point multiply of two
> > single precision numbers, at least when it was being summed into a
double
> >precision variable.  He then explains that IBM changed this in IBSYS 13,
> > anticipating S/360.

> They didn't change the way the 7090 and 7094 computed the result, which
was
> that the LSB of the 36x36 multiply was physically in a register --
> available.  If the programmer had specified that the result of
> single*single was to be placed in a double, then it could be, with no
> further work.  Obviously, single=single*single didn't produce a double
> result.
>
>  > As far as I can tell, the Fortran-66 standard requires a
> > single precision product when multiplying single precision numbers,
> > though I don't know how much impact IBM had on Fortran-66.
>
> I'd appreciate an exact quote on that.  I used to know X3.9-1966 awfully
> well, and I don't recall that the Honeywell gear that I used after the
7094
> and 360 and CDC 3300 had a standards violation in that regard.  My copy is
> filed away, though.

Here is the quote, the last paragraph of section 6.1.  The question is, what
does it mean:

"By use of arithmetic operators other than exponentiation, and admissible
element may be combined with another element of the of the same type, and
the resultant element is of the same type.  Further, an admissible real
element may be combined with a double precision or complex element; the
resultant element is of type double precision or complex, respectively."

>  > The reason why this comment is surprising to me is that S/360
> >(along with S/370 and ESA/390) do not have a single precision
> > multipy instruction with a single precision product.
> > The product is always double precision.  I have checked the code
> > generated by both the G and H compilers, and both, after
>  > generating the multiply instruction generate at least two more
> > instructions to truncate the product to single precision.

> I wonder why they bother.  If the result is to be stored, the extra
> precision goes away.  If it's to continue in a single-precision
> computation, it doesn't help much -- only the first bit of the LSBs has
any
> effect on the next multiplication, and rounding introduces a possible
error
> of the same magnitude.

If the result is stored in a double precision variable, or summed into such
variable, it doesn't.

>  > As far as I know, Fortran still requires floating point multiply to
> >generate a result with no more precision than the arguments that went
> > into it.  W.Kahan describes well the advantage of summing a double
> >precision product of single precision arguments,
> > yet the Fortran standard doesn't seem to follow it.
>  > (Someone with a more recent standard can comment on this.)

> Honeywell had what are called overlength registers, as Tim Prince has
> mentioned.  A 36-bit float was actually an 8-bit exponent and a 28-bit
> fraction.  On loading these into registers, they went into an 8-bit
> E-register, and a 36-bit (machine word-length) register.  A following
> multiplication, for example, produced a normalized result by updating the
> Exponent register, with overlength arithmetic product in the full-length
> register.  That overlength product or other intermediate result was
> maintained so long as the computation did not send a result to storage,
but
> a floating-store rounded the value properly and put it back into a 36-bit
> (8+28 bit) float.  It was a terrific system and produced very good, very
> acceptable to programmers, results.  I'm frankly skeptical that this
> violated X3.9-1966.

>  > I would, in the end, like to thank W.Kahan for his unceasing drive to
>  > improve the results of floating point hardware and software, and hope
> > for only constructive discussion of this problem.

> That cite, so the rest of us could look at the article to which you're
> referring, would be real nice.
>
> I had to respond to government (US national government) queries about the
> Honeywell overlength registers on a regular basis, and I used to include
> with my replies to specific questions a copy of a paper by W. D. Maurer
> (I'm pretty sure that's the one), from Communications of the ACM.  His
> paper was on the results of floating point arithmetic operations on many
> platforms, with the overlength issue a minor mention near the end of the
> paper.  I don't have a precise cite but the paper came out in roughly
1972.
>  It may just have been my magic personality or my renowned writing skills,
> but I like to think that we pretty well settled all the concerns they had
> about Honeywell's somewhat unusual arithmetic arrangement.  We were never
> disqualified from a bid, and I don't recall any case in which we were
> either marked-down or questioned further, once the potential customer had
a
> chance to review our specific comments and the paper from CACM.

-- glen

"Glen Herrmannsfeldt" <gah@ugcs.caltech.edu> writes:

> As far as I know, Fortran still requires floating point multiply to generate
> a result with no more precision than the arguments that went into it.

I disagree.  I don't think that any version of the standard makes such a
requirement.  The standard does not specify accuracy requirements.  In
my opinion, you are getting into the question of computational versus
mathematical differences.  I just got up and don't feel like dragging out the
latest standard to find the relevant citation (much less finding the
similar citation in other versions), but I'm sure many here will know
the one I'm referring to.

The fact that the standard describes  (a*b)+(c*d) interms of single
precision (when a,b,c,d are all single) does not mean that the result of
that computation has to be bit-by-bit identical with a naive implementation
of that expression.

If you go down that path, before long you'll probably conclude that the
standard disallows guard bits.  I simply do not believe this to be true
(or ever to have been true).

As someone else pointed out, the standard (as of f77, I think it was)
provides the DPROD intrinsic specifically for the purpose of doing a
double precision product of single precision operands.  Let me note
on a related matter that the standard does *NOT* require the result of
dprod(a,b) to actually be any more accurate than a*b.  It has to be
represented in double precision, but it doesn't have to meet any particular
accuracy requirements.  Of course, a user is likely to expect it to,
but I think that's one of thoes quality of implementation issues.  If
the vendor doesn't feel like bothering to do a good job of it, the
standard isn't the tool that will address that problem in this case.

-- 
Richard Maine
email: my last name at domain
domain: isomedia dot com

On Sat, 27 Sep 2003 00:17:13 GMT, "Glen Herrmannsfeldt"
<gah@ugcs.caltech.edu> wrote:

<snip>
>
>PL/I has the MULTIPLY function where one can specify the desired result
>precision.  At least one compiler I know of expands both operands before
>doing a double precision multiply when double precision is desired.
>
>If someone wanted to comment about other hardware or compilers in use in the
>1960's and 1970's, that might also be useful.
>

When I took numerical methods as a graduate student (mid-seventies), I
knew only PL/I and 360 Assembler.  Everybody was running on IBM 360
(or whatever the model was at the time).  My classmates, who all used
Fortran IV (while I was coding in PL/I), invariably got different
answers (in the least significant digits) doing the same problem on
the same machine.  Since I was not majoring in computer science,
neither the professor nor I had the competence or the interest to
explore the details, although I might take a different attitude today.

I got my thesis written without ever learning Fortran, but my
graduation into the world of employment put an end to my dogged
loyalty to PL/I, a language I still miss.

I'd like to propose the eccentric view that having little oddities of
floating point arithmetic floating around isn't necessarily a bad
thing.  A rigidly enforced monoculture of floating point arithmetic
could mask problems that otherwise might surface in a slightly less
precise world.

RM

"Robert Myers" <rmyers@rustuck.com> wrote in message
news:iu6fnvch0g09q9f335tt89ogfcj5jm2ovh@4ax.com...
> On Sat, 27 Sep 2003 00:17:13 GMT, "Glen Herrmannsfeldt"
> <gah@ugcs.caltech.edu> wrote:
>
> <snip>
> >
> >PL/I has the MULTIPLY function where one can specify the desired result
> >precision.  At least one compiler I know of expands both operands before
> >doing a double precision multiply when double precision is desired.
> >
> >If someone wanted to comment about other hardware or compilers in use in
the
> >1960's and 1970's, that might also be useful.
> >
>
> When I took numerical methods as a graduate student (mid-seventies), I
> knew only PL/I and 360 Assembler.  Everybody was running on IBM 360
> (or whatever the model was at the time).  My classmates, who all used
> Fortran IV (while I was coding in PL/I), invariably got different
> answers (in the least significant digits) doing the same problem on
> the same machine.  Since I was not majoring in computer science,
> neither the professor nor I had the competence or the interest to
> explore the details, although I might take a different attitude today.

Since my post on Fortran G and H, I have also run the same test on PL/I (F).
That might be the compiler you used.

For   C=C+A*B; single precision A and B, and double precision C, it does the
same thing Fortran did.

It would be interesting to see the results from newer Fortran and PL/I
compilers.

-- glen

In article <kmRdb.626671$YN5.452649@sccrnsc01>,
"Glen Herrmannsfeldt" <gah@ugcs.caltech.edu> writes:
|> 
|> Since my post on Fortran G and H, I have also run the same test on PL/I (F).
|> That might be the compiler you used.
|> 
|> For   C=C+A*B; single precision A and B, and double precision C, it does the
|> same thing Fortran did.
|> 
|> It would be interesting to see the results from newer Fortran and PL/I
|> compilers.

Yes, but you haven't posted WHAT it did!  Given that it doesn't need
to do anything special, I am puzzled over what instructions it inserts.
Can you post the generated code?


Regards,
Nick Maclaren.

"Nick Maclaren" <nmm1@cus.cam.ac.uk> wrote in message
news:bl8q1g$iv7$1@pegasus.csx.cam.ac.uk...
>
> In article <kmRdb.626671$YN5.452649@sccrnsc01>,
> "Glen Herrmannsfeldt" <gah@ugcs.caltech.edu> writes:
> |>
> |> Since my post on Fortran G and H, I have also run the same test on PL/I
(F).
> |> That might be the compiler you used.
> |>
> |> For   C=C+A*B; single precision A and B, and double precision C, it
does the
> |> same thing Fortran did.
> |>
> |> It would be interesting to see the results from newer Fortran and PL/I
> |> compilers.
>
> Yes, but you haven't posted WHAT it did!  Given that it doesn't need
> to do anything special, I am puzzled over what instructions it inserts.
> Can you post the generated code?

I thought it was somewhere in the thread, but maybe not.

For  C=C+A*B, single precision A, B, and double precision C:

   LE 0,A   load single precision A
   ME 0,B   Multiply by B generating double product in float register 0
   SDR 2,2   Clear double float register 2
   LER 2,0    copy single value from reg 0 to reg 2
   AD 2,C     add double C
   STD 2,C    store double into C

The SDR and LER are used to clear the low half of the product.

-- glen

| Richard Maine wrote:
|> Ron Shepard wrote:
|> Are the exponents converted automatically with these instructions?
|> If they aren't, then it seems odd to be able to convert between a
|> double precision floating point value to a single precision value
|> just by loading the value into a register.
|
| Add usual caveat about my memory, but the exponents don't *NEED*
| converting.  I recall thinking it sort of strange the first time I saw
| floating point formats where the exponents weren't the same for single
| and double.  I've subsequently gotten used to it.  Sounds to me like
| you are so used to it that you didn't think of the other possibility.  :-)

What if the number isn't normalized ? ___________________________Gerard S.

"GerardS" <GerardS@PrairieTech.Net> writes:

> | Richard Maine wrote:
> |> Ron Shepard wrote:
> |> Are the exponents converted automatically with these instructions?
> |> If they aren't, then it seems odd to be able to convert between a
> |> double precision floating point value to a single precision value
> |> just by loading the value into a register.
> |
> | Add usual caveat about my memory, but the exponents don't *NEED*
> | converting.  I recall thinking it sort of strange the first time I saw
> | floating point formats where the exponents weren't the same for single
> | and double.  I've subsequently gotten used to it.  Sounds to me like
> | you are so used to it that you didn't think of the other possibility.  :-)
> 
> What if the number isn't normalized ? ___________________________Gerard S.

Why would it make any difference?  It doesn't that I know of.  Do you
have something in specific in mind?

Yes, numbers that are not normalized are prone to lose precision if
you don't normalize them before operations.  In particular, if you
convert a nonnormalized double to a single, you'll presumably get
a nonnormalized single, which will have worse precision than the
horrid IBM single format does anyway.

If one goes around doing operations on nonnormalized data, you are
likely to lose precision in general.  Nothing unique to
single-double conversion here that I see.

Again, I ackowledge that my memory may be confusing systems.  I
certainly recall (with much pain) the IBM single precision floating
point.  I get a little hazy about the details of the machine
instruction set, though.  I did once do some BAL coding - even took a
formal course in it, but it's been a long time (nearly 30 years) and I
never did a lot of it.

-- 
Richard Maine                       |  Good judgment comes from experience;
email: my first.last at org.domain  |  experience comes from bad judgment.
org: nasa, domain: gov              |        -- Mark Twain

In article <uebrt3frgh.fsf@altair.dfrc.nasa.gov>,
Richard Maine  <nospam@see.signature> wrote:
>"GerardS" <GerardS@PrairieTech.Net> writes:
>> | Richard Maine wrote:
>> |> Ron Shepard wrote:
>>
>> |> Are the exponents converted automatically with these instructions?
>> |> If they aren't, then it seems odd to be able to convert between a
>> |> double precision floating point value to a single precision value
>> |> just by loading the value into a register.
>> |
>> | Add usual caveat about my memory, but the exponents don't *NEED*
>> | converting.  I recall thinking it sort of strange the first time I saw
>> | floating point formats where the exponents weren't the same for single
>> | and double.  I've subsequently gotten used to it.  Sounds to me like
>> | you are so used to it that you didn't think of the other possibility.  :-)
>> 
>> What if the number isn't normalized ? ___________________________Gerard S.
>
>Why would it make any difference?  It doesn't that I know of.  Do you
>have something in specific in mind?

The only issue is that unusual amounts of precision can be lost
(including converting a normal range non-zero to a zero), but that
is not really a problem.

>Yes, numbers that are not normalized are prone to lose precision if
>you don't normalize them before operations.  In particular, if you
>convert a nonnormalized double to a single, you'll presumably get
>a nonnormalized single, which will have worse precision than the
>horrid IBM single format does anyway.

Yes.  If I recall, most operations prenormalised, but division didn't,
so using unnormalised numbers in divide was BAD NEWS.

>If one goes around doing operations on nonnormalized data, you are
>likely to lose precision in general.  Nothing unique to
>single-double conversion here that I see.

I am 99% certain that there was nothing special there.

>Again, I ackowledge that my memory may be confusing systems.  I
>certainly recall (with much pain) the IBM single precision floating
>point.  I get a little hazy about the details of the machine
>instruction set, though.  I did once do some BAL coding - even took a
>formal course in it, but it's been a long time (nearly 30 years) and I
>never did a lot of it.

It isn't.  My memory is fading, but it is still good enough to be
certain of that, and I have a POP in my office.


Regards,
Nick Maclaren.

On Mon, 29 Sep 2003 11:51:52 -0500 in alt.folklore.computers,
"GerardS" <GerardS@PrairieTech.Net> wrote:

>| Richard Maine wrote:
>|> Ron Shepard wrote:
>|> Are the exponents converted automatically with these instructions?
>|> If they aren't, then it seems odd to be able to convert between a
>|> double precision floating point value to a single precision value
>|> just by loading the value into a register.
>|
>| Add usual caveat about my memory, but the exponents don't *NEED*
>| converting.  I recall thinking it sort of strange the first time I saw
>| floating point formats where the exponents weren't the same for single
>| and double.  I've subsequently gotten used to it.  Sounds to me like
>| you are so used to it that you didn't think of the other possibility.  :-)
>
>What if the number isn't normalized ? ___________________________Gerard S.

IBM FP does not support denormals. 

Thanks. Take care, Brian Inglis 	Calgary, Alberta, Canada
-- 
Brian.Inglis@CSi.com 	(Brian dot Inglis at SystematicSw dot ab dot ca)
    fake address		use address above to reply

In article <cevgnvkev8rmbb4q6qtmnp1mm6a26ee17n@4ax.com>,
Brian Inglis  <Brian.Inglis@SystematicSw.ab.ca> wrote:
>>
>>What if the number isn't normalized ? ___________________________Gerard S.
>
>IBM FP does not support denormals. 

So?  It most definitely supports unnormalised numbers, which are very
different from IEEE denormals, but may be introduced into IEEE 754R
(for decimal only).


Regards,
Nick Maclaren.

Brian Inglis <Brian.Inglis@SystematicSw.ab.ca> writes:

> On Mon, 29 Sep 2003 11:51:52 -0500 in alt.folklore.computers,
> "GerardS" <GerardS@PrairieTech.Net> wrote:

> >What if the number isn't normalized ? ___________________________Gerard S.
> 
> IBM FP does not support denormals. 

Denormals are not the same thing as unnormalized numbers.  Not the
same thing at all.  Correct - the IBM FP format has no such thing as
denormals, but it does have unnormalized numbers.  Conversely, the
IEEE FP has denormals, but does not support storage of unnormalized
numbers.  The relationship between the two is tenuous.  IEEE
denormals are for very small values, almost underflowing.  IBM
unnormalized numbers can happen for any magnitude.

The IBM FP format does support unnormalized numbers... though the
assembly-language programmer is well advised to be careful about them.

-- 
Richard Maine                       |  Good judgment comes from experience;
email: my first.last at org.domain  |  experience comes from bad judgment.
org: nasa, domain: gov              |        -- Mark Twain

"Nick Maclaren" <nmm1@cus.cam.ac.uk> wrote in message
news:bl9ild$abv$1@pegasus.csx.cam.ac.uk...
>
> In article <xrXdb.631209$Ho3.121609@sccrnsc03>,
> "Glen Herrmannsfeldt" <gah@ugcs.caltech.edu> writes:
> |>
> |> I thought it was somewhere in the thread, but maybe not.
> |>
> |> For  C=C+A*B, single precision A, B, and double precision C:
> |>
> |>    LE 0,A   load single precision A
> |>    ME 0,B   Multiply by B generating double product in float register 0
> |>    SDR 2,2   Clear double float register 2
> |>    LER 2,0    copy single value from reg 0 to reg 2
> |>    AD 2,C     add double C
> |>    STD 2,C    store double into C
> |>
> |> The SDR and LER are used to clear the low half of the product.
>
> If my memory is right, those two are NOT done as part of the
> multiplication, but are the conversion from single to double
> precision.  I.e. you would also get them in:
>
>     C=C+A, single precision A, and double precision C.

Well, you would get SDR 2,2 and LE 2,A in that case.  The compiler does have
to keep track of which registers have single and which have double values.
Leaving low bits from previous calculations in would seem a bad idea.
Though according to other posts it may be legal in Fortran.

-- glen

In article <Rh0eb.632709$Ho3.121617@sccrnsc03>,
Glen Herrmannsfeldt <gah@ugcs.caltech.edu> wrote:
>
>Well, you would get SDR 2,2 and LE 2,A in that case.  The compiler does have
>to keep track of which registers have single and which have double values.
>Leaving low bits from previous calculations in would seem a bad idea.
>Though according to other posts it may be legal in Fortran.

It certainly is legal, provided they are from an appropriate calculation,
but that is not what is going on here.

The specification of IBM floating-point register operations is precisely
that the unused low-order bits are ignored.  Those bits are therefore
NOT being left in the value; they are being left OUTSIDE the value.
That is why it is quite safe for single precision multiply to create
them and to leave them where they are.

Similarly, conversion from single to double precision has to clear those
bits, whether the operand being converted has been created by multiply
or anything else.


Regards,
Nick Maclaren.

"GerardS" <GerardS@PrairieTech.Net> wrote in message
news:vngokjk8nourca@corp.supernews.com...
> | Richard Maine wrote:
> |> Ron Shepard wrote:
> |> Are the exponents converted automatically with these instructions?
> |> If they aren't, then it seems odd to be able to convert between a
> |> double precision floating point value to a single precision value
> |> just by loading the value into a register.
> |
> | Add usual caveat about my memory, but the exponents don't *NEED*
> | converting.  I recall thinking it sort of strange the first time I saw
> | floating point formats where the exponents weren't the same for single
> | and double.  I've subsequently gotten used to it.  Sounds to me like
> | you are so used to it that you didn't think of the other possibility.
:-)
>
> What if the number isn't normalized ? ___________________________Gerard S.

For hex floating point, the exponent is 7 bits for short, long, and extended
formats, allowing 16**-64 up to 16**63, or something like 10**-78 up to
10**75.

Unnormalized values are also allowed in mathematical operations.  For
addition and subtraction, the prenormalization (aligning the hexadecimal
point) is done based on the exponent supplied.  The Fortran AINT function
consists of adding 0*16**7, and DINT of adding 0*16**15.  Prenormalization
will shift out fractional digits, zero will be added, and post normalization
will shift it back to properly normalized form.

Multipy, divide, and square root prenormalize the the operand if it isn't
normalized.  (This is from ESA/390.)

http://publibfp.boulder.ibm.com/cgi-bin/bookmgr/BOOKS/dz9ar006/18.1.2

-- glen

"Nick Maclaren" <nmm1@cus.cam.ac.uk> wrote in message
news:bla4k2$nas$1@pegasus.csx.cam.ac.uk...
> In article <Rh0eb.632709$Ho3.121617@sccrnsc03>,
> Glen Herrmannsfeldt <gah@ugcs.caltech.edu> wrote:
> >
> >Well, you would get SDR 2,2 and LE 2,A in that case.  The compiler does
have
> >to keep track of which registers have single and which have double
values.
> >Leaving low bits from previous calculations in would seem a bad idea.
> >Though according to other posts it may be legal in Fortran.
>
> It certainly is legal, provided they are from an appropriate calculation,
> but that is not what is going on here.
>
> The specification of IBM floating-point register operations is precisely
> that the unused low-order bits are ignored.  Those bits are therefore
> NOT being left in the value; they are being left OUTSIDE the value.
> That is why it is quite safe for single precision multiply to create
> them and to leave them where they are.
>
> Similarly, conversion from single to double precision has to clear those
> bits, whether the operand being converted has been created by multiply
> or anything else.

What I was saying is that if conversion from single to double did not clear
those bits you would get whatever was left there, possibly from an unrelated
calculation.   I believe that would be a bad idea, but some posts by others
seem to indicate that it would be legal Fortran.  Someone else posted:

|Item (6) listed in Section 1.4 of the Fortran standard explicitly
|states that the standard does not specify
|
|    (6)  The physical properties of the representation of quantities and
|         the method of rounding, approximating, or computing numeric
|         values on a particular processor.

The belief is that this statement could include leaving garbage bits in for
low bits.

-- glen

Glen Herrmannsfeldt wrote:
...
> What I was saying is that if conversion from single to double did not clear
> those bits you would get whatever was left there, possibly from an unrelated
> calculation.   I believe that would be a bad idea, but some posts by others
> seem to indicate that it would be legal Fortran.  Someone else posted:
>
>> Item (6) listed in Section 1.4 of the Fortran standard explicitly
>> states that the standard does not specify
>>
>>    (6)  The physical properties of the representation of quantities and
>>         the method of rounding, approximating, or computing numeric
>>         values on a particular processor.
>
> The belief is that this statement could include leaving garbage bits in for
> low bits.

That's probably a correct assessment.  In fact, on the system in question,
where "rounding" results to single precision is done by truncation, it
could be argued that leaving "garbage bits" simulates the missing data!
Yes, I have heard appalling reasoning of that kind.

-- 
J. Giles

"James Giles" <jamesgiles@worldnet.att.net> wrote in message
news:r%0eb.163529$3o3.11771639@bgtnsc05-news.ops.worldnet.att.net...
> Glen Herrmannsfeldt wrote:
> ...
> > What I was saying is that if conversion from single to double did not
clear
> > those bits you would get whatever was left there, possibly from an
unrelated
> > calculation.   I believe that would be a bad idea, but some posts by
others
> > seem to indicate that it would be legal Fortran.  Someone else posted:
> >
> >> Item (6) listed in Section 1.4 of the Fortran standard explicitly
> >> states that the standard does not specify
> >>
> >>    (6)  The physical properties of the representation of quantities and
> >>         the method of rounding, approximating, or computing numeric
> >>         values on a particular processor.
> >
> > The belief is that this statement could include leaving garbage bits in
for
> > low bits.
>
> That's probably a correct assessment.  In fact, on the system in question,
> where "rounding" results to single precision is done by truncation, it
> could be argued that leaving "garbage bits" simulates the missing data!
> Yes, I have heard appalling reasoning of that kind.

Well, since truncating is systematic at puting in zeros, random bits would,
on the average, be closer.

The process of digitizing analog data, such as converting analog audio into
digital form, works best with a process called dithering.  Random noise of
approximately 1 LSB in amplitude is added before sampling.   The left over
bits in registers might not be very random, though.   Using true random
bits, one might even be able to sum into a single precision register and get
good results!

-- glen

Glen Herrmannsfeldt wrote:
> "James Giles" <jamesgiles@worldnet.att.net> wrote in message
> news:r%0eb.163529$3o3.11771639@bgtnsc05-news.ops.worldnet.att.net...
...
>> That's probably a correct assessment.  In fact, on the system in question,
>> where "rounding" results to single precision is done by truncation, it
>> could be argued that leaving "garbage bits" simulates the missing data!
>> Yes, I have heard appalling reasoning of that kind.
>
> Well, since truncating is systematic at puting in zeros, random bits would,
> on the average, be closer.
>
> The process of digitizing analog data, such as converting analog audio into
> digital form, works best with a process called dithering.  Random noise of
> approximately 1 LSB in amplitude is added before sampling.   The left over
> bits in registers might not be very random, though.   Using true random
> bits, one might even be able to sum into a single precision register and get
> good results!

Or not.  The problem is that using "garbage bits" doesn't permit any
useful or reliable analysis of the algorithms.  Clearing the bits does.
Further, clearing the bits allows algorithmic work-arounds (like
Kahan's summation).  Keeping "garbage bits" simply makes that
inoperable.

-- 
J. Giles

On 29 Sep 2003 18:59:14 GMT in alt.folklore.computers,
nmm1@cus.cam.ac.uk (Nick Maclaren) wrote:

>In article <cevgnvkev8rmbb4q6qtmnp1mm6a26ee17n@4ax.com>,
>Brian Inglis  <Brian.Inglis@SystematicSw.ab.ca> wrote:
>>>
>>>What if the number isn't normalized ? ___________________________Gerard S.
>>
>>IBM FP does not support denormals. 
>
>So?  It most definitely supports unnormalised numbers, which are very
>different from IEEE denormals, but may be introduced into IEEE 754R
>(for decimal only).

I ass-u-me-d the poster was asking from an IEEE viewpoint, where
denormals are the only unnormalized form; I also ass-u-me-d, as
you pointed out in another reply, it's pretty obvious there's not
really anything special about possibly losing bits if you operate
on unnormalized numbers. 

Thanks. Take care, Brian Inglis 	Calgary, Alberta, Canada
-- 
Brian.Inglis@CSi.com 	(Brian dot Inglis at SystematicSw dot ab dot ca)
    fake address		use address above to reply

In article <TD0eb.630592$YN5.454239@sccrnsc01>,
"Glen Herrmannsfeldt" <gah@ugcs.caltech.edu> writes:
|> 
|> Unnormalized values are also allowed in mathematical operations.  For
|> addition and subtraction, the prenormalization (aligning the hexadecimal
|> point) is done based on the exponent supplied.  The Fortran AINT function
|> consists of adding 0*16**7, and DINT of adding 0*16**15.  Prenormalization
|> will shift out fractional digits, zero will be added, and post normalization
|> will shift it back to properly normalized form.

That was model dependent.  If I recall, neither IBM's nor Hitachi's
vector units accepted unnormalised operands.

|> Multipy, divide, and square root prenormalize the the operand if it isn't
|> normalized.  (This is from ESA/390.)

That was model dependent, too, if I recall _ I think that multiply
and divide were prone to not prenormalising, at least on some of
the System/360 systems.


Regards,
Nick Maclaren.

In article <Zs2eb.157980$0v4.11847973@bgtnsc04-news.ops.worldnet.att.net>,
"James Giles" <jamesgiles@worldnet.att.net> writes:
|> Glen Herrmannsfeldt wrote:
|> > "James Giles" <jamesgiles@worldnet.att.net> wrote in message
|> > news:r%0eb.163529$3o3.11771639@bgtnsc05-news.ops.worldnet.att.net...
|> ...
|> >> That's probably a correct assessment.  In fact, on the system in question,
|> >> where "rounding" results to single precision is done by truncation, it
|> >> could be argued that leaving "garbage bits" simulates the missing data!
|> >> Yes, I have heard appalling reasoning of that kind.
|> >
|> > Well, since truncating is systematic at puting in zeros, random bits would,
|> > on the average, be closer.
|> >
|> > The process of digitizing analog data, such as converting analog audio into
|> > digital form, works best with a process called dithering.  Random noise of
|> > approximately 1 LSB in amplitude is added before sampling.   The left over
|> > bits in registers might not be very random, though.   Using true random
|> > bits, one might even be able to sum into a single precision register and get
|> > good results!
|> 
|> Or not.  The problem is that using "garbage bits" doesn't permit any
|> useful or reliable analysis of the algorithms.  Clearing the bits does.
|> Further, clearing the bits allows algorithmic work-arounds (like
|> Kahan's summation).  Keeping "garbage bits" simply makes that
|> inoperable.

Yes, definitely.  Despite Kahan's analysis, there are very good
reasons to want true probabilistic rounding.  Several of us have
analysed that it has BETTER properties than nearest rounding (of
any type) in some important operations, like vector summation.
My limited experiments indicated that these benefits mean that it
delivers comparable accuracy to nearest rounding on at least one
of the analyses where the latter does well.

But note that means TRUE probabilistic rounding - accept no
imitations!


Regards,
Nick Maclaren.

"Nick Maclaren" <nmm1@cus.cam.ac.uk> wrote in message
news:blbcjf$lrn$1@pegasus.csx.cam.ac.uk...
>
> In article <Zs2eb.157980$0v4.11847973@bgtnsc04-news.ops.worldnet.att.net>,
> "James Giles" <jamesgiles@worldnet.att.net> writes:
> |> Glen Herrmannsfeldt wrote:
> |> > "James Giles" <jamesgiles@worldnet.att.net> wrote in message
> |> > news:r%0eb.163529$3o3.11771639@bgtnsc05-news.ops.worldnet.att.net...

(snip)

> |> > Well, since truncating is systematic at puting in zeros, random bits
would,
> |> > on the average, be closer.
> |> >
> |> > The process of digitizing analog data, such as converting analog
audio into
> |> > digital form, works best with a process called dithering.  Random
noise of
> |> > approximately 1 LSB in amplitude is added before sampling.   The left
over
> |> > bits in registers might not be very random, though.   Using true
random
> |> > bits, one might even be able to sum into a single precision register
and get
> |> > good results!
> |>
> |> Or not.  The problem is that using "garbage bits" doesn't permit any
> |> useful or reliable analysis of the algorithms.  Clearing the bits does.
> |> Further, clearing the bits allows algorithmic work-arounds (like
> |> Kahan's summation).  Keeping "garbage bits" simply makes that
> |> inoperable.
>
> Yes, definitely.  Despite Kahan's analysis, there are very good
> reasons to want true probabilistic rounding.  Several of us have
> analysed that it has BETTER properties than nearest rounding (of
> any type) in some important operations, like vector summation.
> My limited experiments indicated that these benefits mean that it
> delivers comparable accuracy to nearest rounding on at least one
> of the analyses where the latter does well.
>
> But note that means TRUE probabilistic rounding - accept no
> imitations!

Not so long ago, someone came up with this statement from the standard:

>> Item (6) listed in Section 1.4 of the Fortran standard explicitly
>> states that the standard does not specify
>>
>>    (6)  The physical properties of the representation of quantities and
>>         the method of rounding, approximating, or computing numeric
>>         values on a particular processor.

I had tried to argue that at least the rounding would be expected to be
consistent (repeatable), but the discussion seems to say that it doesn't
need to be.   But yes, I agree that only true probabilistic rounding should
be used, if any is.

It seems that some processors now have a hardware random number generator
using an electronic noise source.  I don't know how fast it is, though.   It
would be interesting to see that work.

-- glen

In article <QEbeb.377802$2x.116648@rwcrnsc52.ops.asp.att.net>,
"Glen Herrmannsfeldt" <gah@ugcs.caltech.edu> writes:
|> 
|> > But note that means TRUE probabilistic rounding - accept no
|> > imitations!
|> 
|> Not so long ago, someone came up with this statement from the standard:
|> 
|> >> Item (6) listed in Section 1.4 of the Fortran standard explicitly
|> >> states that the standard does not specify
|> >>
|> >>    (6)  The physical properties of the representation of quantities and
|> >>         the method of rounding, approximating, or computing numeric
|> >>         values on a particular processor.
|> 
|> I had tried to argue that at least the rounding would be expected to be
|> consistent (repeatable), but the discussion seems to say that it doesn't
|> need to be.   But yes, I agree that only true probabilistic rounding should
|> be used, if any is.

That is my understanding, too.

|> It seems that some processors now have a hardware random number generator
|> using an electronic noise source.  I don't know how fast it is, though.   It
|> would be interesting to see that work.

It would be very interesting.  Generator efficiency isn't a problem,
because you can have (say) a decent shift register generator running
in parallel (e.g. all the time!), and XOR the hardware generator in
with the results of that when you need them.  It wouldn't matter if
the hardware generator updated a word only every 100 cycles, given
a suitably decent pseudorandom generator.

I am not enough of an expert on the details of the methods used in
hardware to be quite certain that there need be NO performance loss,
but I am 99% sure that it would be negligible at worst and there
might even be a slight gain.  Yes, really.

The reason I say the last is that there is no need to calculate a
multiplication to full accuracy - all you need is to generate
enough bits for a decision.  Of course, this means that operations
can take variable times, which is Not Approved Of.  But it might
help, especially for extended precision multiplication.


Regards,
Nick Maclaren.

"Nick Maclaren" <nmm1@cus.cam.ac.uk> wrote in message
news:blbmni$g2$1@pegasus.csx.cam.ac.uk...
>
> In article <QEbeb.377802$2x.116648@rwcrnsc52.ops.asp.att.net>,
> "Glen Herrmannsfeldt" <gah@ugcs.caltech.edu> writes:
> |>
> |> > But note that means TRUE probabilistic rounding - accept no
> |> > imitations!
> |>
[...]

See
DS Parker et al., "Monte Carlo Arithmetic; exploiting randomness in
floating-point arithmetic," Computing in Science and Engineering pp.
58-68 July 2000.

-- 
Ciao,
Gerry T.
______
"Competent engineers rightly distrust all numerical computations and
seek corroboration from alternative numerical methods, from scale
models, from prototypes, from experience... ."  W.V. Kahan.

Brian Inglis <Brian.Inglis@SystematicSw.ab.ca> writes:

> On 29 Sep 2003 18:59:14 GMT in alt.folklore.computers,
> nmm1@cus.cam.ac.uk (Nick Maclaren) wrote:
> 
> >In article <cevgnvkev8rmbb4q6qtmnp1mm6a26ee17n@4ax.com>,
> >Brian Inglis  <Brian.Inglis@SystematicSw.ab.ca> wrote:
> >>>
> >>>What if the number isn't normalized ? ___________________________Gerard S.
> >>
> >>IBM FP does not support denormals. 
> >
> >So?  It most definitely supports unnormalised numbers, which are very
> >different from IEEE denormals, but may be introduced into IEEE 754R
> >(for decimal only).
> 
> I ass-u-me-d the poster was asking from an IEEE viewpoint, where
> denormals are the only unnormalized form;

Oh.  Well, the context was IBM mainframe instructions, which don't
fit that assumption.  At least that's what I thought the context was
unless I got lost.

-- 
Richard Maine
email: my last name at domain
domain: isomedia dot com

On 30 Sep 2003 10:42:58 GMT, nmm1@cus.cam.ac.uk (Nick Maclaren) wrote:

[probabilistic floating point rounding]

>The reason I say the last is that there is no need to calculate a
>multiplication to full accuracy - all you need is to generate
>enough bits for a decision.  Of course, this means that operations
>can take variable times, which is Not Approved Of.  But it might
>help, especially for extended precision multiplication.

I'm no expert on numerical analysis, so tell me if I'm missing
something obvious, but:

Let's suppose you have N bits of precision (N ~= 53 for IEEE double
precision IIRC).

In the base case where rounding truncates; you have N bits of
precision.

If rounding is to nearest, .5 rounding away from zero (the simplest
type of rounding, it requires the consideration of only one extra
bit), you effectively have something like N+0.5 bits of precision in
calculations.

With more complicated rounding, such as to even (IIRC, IEEE supports
modes like this, at considerable complexity penalty), or
probabilistic, you effectively have something like N+0.5+epsilon bits.

Now my question is: why is it worth so much trouble to get an extra
epsilon of a bit, when by using e.g. Intel's 80 bit format for
temporaries, you get something like an extra 16 bits?

-- 
"Sore wa himitsu desu."
To reply by email, remove
the small snack from address.
http://www.esatclear.ie/~rwallace

In article <3f7a622d.113216013@news.eircom.net>,
wallacethinmintr@eircom.net (Russell Wallace) writes:
|> 
|> >The reason I say the last is that there is no need to calculate a
|> >multiplication to full accuracy - all you need is to generate
|> >enough bits for a decision.  Of course, this means that operations
|> >can take variable times, which is Not Approved Of.  But it might
|> >help, especially for extended precision multiplication.
|> 
|> I'm no expert on numerical analysis, so tell me if I'm missing
|> something obvious, but:
|> 
|> Let's suppose you have N bits of precision (N ~= 53 for IEEE double
|> precision IIRC).
|> 
|> In the base case where rounding truncates; you have N bits of
|> precision.

OK.

|> If rounding is to nearest, .5 rounding away from zero (the simplest
|> type of rounding, it requires the consideration of only one extra
|> bit), you effectively have something like N+0.5 bits of precision in
|> calculations.

Well, maybe.

|> With more complicated rounding, such as to even (IIRC, IEEE supports
|> modes like this, at considerable complexity penalty), or
|> probabilistic, you effectively have something like N+0.5+epsilon bits.

Not really, but explaining why takes a while.  It isn't important,
anyway.

|> Now my question is: why is it worth so much trouble to get an extra
|> epsilon of a bit, when by using e.g. Intel's 80 bit format for
|> temporaries, you get something like an extra 16 bits?

Because of the order of error build up.  Consider the simple summation
of a vector of length L:

Truncation is almost always biassed, so that the relative error will
be of order L*eps/2.  Most roundings are not much biassed for well
mixed inputs and so the error will be of order sqrt(L)*eps/4.  A
major gain.

In the case where the initial value of the sum is A and you are
adding small numbers very near B, where the variation in B is less
than A*eps/L, ANY true rounding will have error of order L*eps/4
(though with a much smaller factor than for truncation).  True
probabilistic rounding, however, maintains an order sqrt(L)*eps/2.8.

Now, if you use extra precision for the sum, this ceases to be a
problem, but that applies only when you know WHERE in your program
the accuracy loss will occur and you use the extended precision
appropriately.  If you use extra precision, then you are merely
reducing eps - sometimes this will eliminate the problem, sometimes
it won't.

This is why hardware designers go to great trouble (and it IS great
trouble) to support rounding.  IBM originally omitted from the
System/360 because doing so can save a complete cycle on floating
point operations.  Note that the problem with rounding is twofold:

    1) You have to maintain at least one extra bit (two for IEEE)
to know how to round.

    2) You may have to renormalise or reraise an exception after
the rounding (and this is where the extra cycle goes, in general).

What I and a few other people have said over many years (and I am
a mere follower) is that the extra cost of PERMITTING a true
probabilistic rounding mode is negligible, and there are certain
calculations for which it is vastly better.  Apparently it is
starting to be adopted in graphics.


Regards,
Nick Maclaren.

"Russell Wallace" <wallacethinmintr@eircom.net> wrote in message
news:3f7a622d.113216013@news.eircom.net...
> On 30 Sep 2003 10:42:58 GMT, nmm1@cus.cam.ac.uk (Nick Maclaren) wrote:
>
> [probabilistic floating point rounding]
>
> >The reason I say the last is that there is no need to calculate a
> >multiplication to full accuracy - all you need is to generate
> >enough bits for a decision.  Of course, this means that operations
> >can take variable times, which is Not Approved Of.  But it might
> >help, especially for extended precision multiplication.
>
> I'm no expert on numerical analysis, so tell me if I'm missing
> something obvious, but:
>
> Let's suppose you have N bits of precision (N ~= 53 for IEEE double
> precision IIRC).
>
> In the base case where rounding truncates; you have N bits of
> precision.
>
> If rounding is to nearest, .5 rounding away from zero (the simplest
> type of rounding, it requires the consideration of only one extra
> bit), you effectively have something like N+0.5 bits of precision in
> calculations.

I am not sure that I understand it so well, but I saw an explanation that to
get proper rounding for division with an N bit quotient, it is necessary to
calculate 2N quotient bits.    I believe that this is related to the
statement about double rounding, where rounding to >=2N bits, and then
subsequently to N bits gets the same result as rounding directly to N bits.
Intermediate rounding to <2N bits before rounding to N bits can give the
wrong answer.

The statement about division is in the IBM Systems Journal article on
extended precision floating point, and the DXR simulation routine.

-- glen

-- glen

> In the case where the initial value of the sum is A and you are
> adding small numbers very near B, where the variation in B is less
> than A*eps/L, ANY true rounding will have error of order L*eps/4
> (though with a much smaller factor than for truncation).  True
> probabilistic rounding, however, maintains an order sqrt(L)*eps/2.8.

This sounds like the same argument for dithering in analog-to-digital
conversion...?

	Jan

In article <3F7AA528.5422AAD9@mediasec.de>, Jan C. =?iso-8859-1?Q?Vorbr=FCggen?= <jvorbrueggen@mediasec.de> writes:
|> > In the case where the initial value of the sum is A and you are
|> > adding small numbers very near B, where the variation in B is less
|> > than A*eps/L, ANY true rounding will have error of order L*eps/4
|> > (though with a much smaller factor than for truncation).  True
|> > probabilistic rounding, however, maintains an order sqrt(L)*eps/2.8.
|> 
|> This sounds like the same argument for dithering in analog-to-digital
|> conversion...?

Yup :-)


Regards,
Nick Maclaren.

In article <HQweb.382949$2x.123651@rwcrnsc52.ops.asp.att.net>,
"Glen Herrmannsfeldt" <gah@ugcs.caltech.edu> writes:
|> 
|> I am not sure that I understand it so well, but I saw an explanation that to
|> get proper rounding for division with an N bit quotient, it is necessary to
|> calculate 2N quotient bits.    I believe that this is related to the
|> statement about double rounding, where rounding to >=2N bits, and then
|> subsequently to N bits gets the same result as rounding directly to N bits.
|> Intermediate rounding to <2N bits before rounding to N bits can give the
|> wrong answer.

Yes, roughly, though it isn't so simple.  But the latter statements
are definitely correct.

|> The statement about division is in the IBM Systems Journal article on
|> extended precision floating point, and the DXR simulation routine.

It predates that, and is well-known.


The rough summary is that addition, subtraction, comparison and so
on need only a fixed number of extra bits to support any sane form
of rounding, where the number is typically 1 or 2 on top of the
guard bits.

Multiplication is a real pain, because it really needs to be done
to a full NxN -> 2N form, and THEN rounded.  There are optimisations,
of course.

Division is definitely tricky, because it depends on how you do it.
I don't know the methods used in hardware, but the various forms of
repeated subtraction mean that you need calculate only N bits, but
need to do that to 2N precision.  Newton-Raphson is similar.

It is fairly easy to prove that the issues are close to fundamental,
and that ANY accurate multiplication or division method will need
the extra precision.  A long time back (1980s?), R.P. Brent did a
PROPER complexity analysis of the basic operations, and showed why
multiplication was more complex than addition.  Division is well
known to be a simple factor of about 2-3 more expensive than
multiplication, of course, for large numbers of bits.

None of this changes when considering probabilistic versus
deterministic rounding.


Regards,
Nick Maclaren.

On 1 Oct 2003 07:38:25 GMT, nmm1@cus.cam.ac.uk (Nick Maclaren) wrote:

[explanation snipped]

>This is why hardware designers go to great trouble (and it IS great
>trouble) to support rounding.

So basically, for an individual operation it indeed doesn't much
matter, but for large calculations you want to use a rounding mode
chosen so that the sum of all the errors will average out as close as
possible to zero. Okay, that makes sense, thanks.

-- 
"Sore wa himitsu desu."
To reply by email, remove
the small snack from address.
http://www.esatclear.ie/~rwallace

> |> This sounds like the same argument for dithering in analog-to-digital
> |> conversion...?
> 
> Yup :-)

Hey, that's good: I've finally understood something about probabilistic
rouding - and from an explanation by you to boot 8-)!!

	Jan

Jan C. Vorbr�ggen wrote:
>>In the case where the initial value of the sum is A and you are
>>adding small numbers very near B, where the variation in B is less
>>than A*eps/L, ANY true rounding will have error of order L*eps/4
>>(though with a much smaller factor than for truncation).  True
>>probabilistic rounding, however, maintains an order sqrt(L)*eps/2.8.
> 
> 
> This sounds like the same argument for dithering in analog-to-digital
> conversion...?

It is _exactly_ analogous to that approach, it avoids biasing all the 
way down to below the actual precision of the sum.

Terje

-- 
- <Terje.Mathisen@hda.hydro.com>
"almost all programming can be viewed as an exercise in caching"

Nick> Apparently it [true probabilistic rounding] is starting to be
Nick> adopted in graphics.

Whoa!

The graphics guys at SGI soured on MADDs when they found out
that fused multiply-add instructions round differently than two
separate operations.  It turns out that you have to calculate the
exact same screen position twice for a vertex which is part of two
different triangles.  The bug they found was that the vertex went
down two different code paths, and the two paths rounded
differently, in this case because one path used multiply and then
add, and the other used a fused multiply-add.

Graphics guys are really into repeatability.  They'll get their
performance by chucking precision or accuracy first.

In article <45022fc8.0310021010.4aaf0e5e@posting.google.com>,
Iain McClatchie <iain-3@truecircuits.com> wrote:
>
>Nick> Apparently it [true probabilistic rounding] is starting to be
>Nick> adopted in graphics.
>
>Whoa!
>
>The graphics guys at SGI soured on MADDs when they found out
>that fused multiply-add instructions round differently than two
>separate operations.  It turns out that you have to calculate the
>exact same screen position twice for a vertex which is part of two
>different triangles.  The bug they found was that the vertex went
>down two different code paths, and the two paths rounded
>differently, in this case because one path used multiply and then
>add, and the other used a fused multiply-add.
>
>Graphics guys are really into repeatability.  They'll get their
>performance by chucking precision or accuracy first.

If they code like that, they will also be sacrificing performance!

The first point is that you don't get repeatability in combination
with significant optimisation - even without that problem.  The
second one is, when calculating a vertex's coordinates and consistency
is critical, do it once and save it.  This can be essential for things
like Dirichlet tesselation, which are seriously numerically unstable.

I was actually referring to the use in dithering, which is somewhere
that true probabilistic rounding is much superior to any deterministic
rounding.

What is probably the case is that anyone who relies on absolutely
consistent floating-point results is living dangerously, and may be
deluding himself.  I know of at least half a dozen ways where the
assumption will fail, even on consistent IEEE-like architectures.
And considerably more on older ones.

That is part of the basis for my assertion that probabilistic rounding
is good because it forces programmers to realise that floating-point 
arithmetic does not give precisely predictable results.


Regards,
Nick Maclaren.

Iain McClatchie wrote:

> Nick> Apparently it [true probabilistic rounding] is starting to be
> Nick> adopted in graphics.
> 
> Whoa!
> 
> The graphics guys at SGI soured on MADDs when they found out
> that fused multiply-add instructions round differently than two
> separate operations.  It turns out that you have to calculate the
> exact same screen position twice for a vertex which is part of two
> different triangles.  The bug they found was that the vertex went
> down two different code paths, and the two paths rounded
> differently, in this case because one path used multiply and then
> add, and the other used a fused multiply-add.
> 
> Graphics guys are really into repeatability.  They'll get their
> performance by chucking precision or accuracy first.

You're absolutely right.

It really doesn't have to be mathematically perfect, but it does have to 
be repeatable.

id Software, in the form of my friend Mike Abrash, spent a few weeks 
catching an intermittent problem of exactly this kind:

Every once in a (long) while, such a duplicated calculation of screen 
coordinates along a vertex would turn up two different pixels, and 
there's be a one-frame-long single-pixel glitch in the screen.

This was only observed a few (3-4?) times, and only on one particular 
machine.

After a _lot_ of code inspections, logging, inserted debugging 
statements etc, they still couldn't come up with a repeatable scenario.

That was when they disassembled the cpu/cooling fan combo, and 
discovered that their vendor had delivered an overclocked (by just 11%, 
90 MHz -> 100 MHz) cpu. After a quick reset of the relevant motherboard 
switches, everything was OK.

So remember: Repeatability isn't an option, it's an absolute requirement.

Terje
-- 
- <Terje.Mathisen@hda.hydro.com>
"almost all programming can be viewed as an exercise in caching"

Nick> ...when calculating a vertex's coordinates and consistency
Nick> is critical, do it once and save it.

They can't.  DRAM pins cost more than multipliers.  As for optimization,
they just have to control their optimizers so that if two different
code paths have to generate identical results, they do.  That's just
part of graphics engineering.

Nick> I was actually referring to the use in dithering, which is somewhere
Nick> that true probabilistic rounding is much superior to any deterministic
Nick> rounding.

Sounds good.  That's really FP->int.

Which is a great opportunity to get off topic.

The human retina has loads more dynamic range than 1:256.  The analog
connection between a video board and the monitor has limited dynamic
range, but now that we've got DVI cables, the limiting factor in dynamic
range is all in one unit: the drivers and liquid crystal or electron guns
and phosphors of the displays.  Surely it's time for more dynamic range.
And since the retina, at least, has nonlinear response, perhaps we need
floating-point DACs driving these displays?

On 2 Oct 2003 18:36:58 GMT, Nick Maclaren <nmm1@cus.cam.ac.uk> wrote:
> What is probably the case is that anyone who relies on absolutely
> consistent floating-point results is living dangerously, and may be
> deluding himself.  I know of at least half a dozen ways where the
> assumption will fail, even on consistent IEEE-like architectures.
> And considerably more on older ones.
> 
> That is part of the basis for my assertion that probabilistic rounding
> is good because it forces programmers to realise that floating-point 
> arithmetic does not give precisely predictable results.

I would like true pseudoprobabilistic rounding, one with the
same statisitics as probabilistic rounding, but determinsitically
repeatable if need be.  Gaussian rounding would be better still.  :) 

In dynamical systems at least there almost always is (in reality)
some kind of tiny probabilistic forces, e.g. thermodynamic fluctuations
on a macroscopic object.

Such rounding in calculations would make numerical results which are
ore like what would "physically happen", in the sense that orbits that
arise are likely to be ones on the "natural" physical measure.

The issue of measure of results is in fact tricky and difficult
in theoretical dynamical systems, but everybody agrees that there
is some kind of "natural" physical measure out there.