On 17 Apr 2024 at 15:08:01 CEST, "PMF" <PMF> wrote:
> In lxf64/ntf64 I use an external C library, specifically fdlibm53.
>
> This claims to be within 1 ulp correct. My testing also suggests this.
> fdlibm looks to be the base for most other math libraries.
> I could have used libm from gcc but I wanted the same code for both
> Linux and Windows.
>
> In lxf/ntf I use the 387 fp stack. I think this was a wrong decision.
> 8 stack items is to low to be useful for anything more complicated.
> complex numbers is an example that quickly eats all 8 stack items.

I share your pain. In 2023, I spent far too long trying to satisfy FP users.
For VFX64, we ended up providing support for
1) 8 item 387 FP, internal FP stack,
2) 387 FP and FP stack in memory,
3) SSE2 FP and FP stack in memory.

Conversion of FP parameters to/from OS call interfaces is done
automagically by the EXTERN: declarations for system calls and
callbacks.

For a significant proportion of our users, 64 bit FP is inadequate,
and ony 64 bit FP is available on CPUs other than AMD64/x64
unless you are using rare and expensive hardware.

In retrospect, if I were doing this again I would standardise on an
external double-double library (about 106 bits). In most cases that
we encounter, the desire for 387 FP is to gain the extra precision.
Since very few CPUs support quad precision natively, the most
obvious solution is a double-double library.

Stephen

--
Stephen Pelc, stephen@vfxforth.com
MicroProcessor Engineering, Ltd. - More Real, Less Time
133 Hill Lane, Southampton SO15 5AF, England
tel: +44 (0)78 0390 3612, +34 649 662 974
http://www.mpeforth.com
MPE website
http://www.vfxforth.com/downloads/VfxCommunity/
downloads

On 4/17/24 08:08, PMF wrote:
> Krishna Myneni wrote:
>
>> On 4/15/24 20:10, minforth wrote:
>>> I would be surprised to get a SIGFPE interrupt from x87 stack ops.
>>> X87 mode has long been deprecated and replaced by SSE2.
>
>> Maybe I'm missing something, but SSE2 does not seem to have anything
>> beyond basic floating point arithmetic e.g. transcendental functions.
>
> Yes that is right. There is a sqrt but nothing more.
>
> In lxf64/ntf64 I use an external C library, specifically fdlibm53.
>
> This claims to be within 1 ulp correct. My testing also suggests this.
> fdlibm looks to be the base for most other math libraries.
> I could have used libm from gcc but I wanted the same code for both
> Linux and Windows.
>
> In lxf/ntf I use the 387 fp stack. I think this was a wrong decision.
> 8 stack items is to low to be useful for anything more complicated.
> complex numbers is an example that quickly eats all 8 stack items.
> Best Regards
> Peter Fälth

Thank you for the clarification. I disassembled the GNU C Math Library
(64-bit libm.so.6) and had a look at the code. It is a huge file. The
long double functions are coded for x87, while functions for float,
double, and float128 use the SSE instructions.

This implies that if x87 is deprecated on 64-bit systems, then 80-bit
floats are also deprecated?

Does fdlibm53 support 80-bit floats and float128?

--
Krishna

On 4/18/24 07:09, Stephen Pelc wrote:
....
> In retrospect, if I were doing this again I would standardise on an
> external double-double library (about 106 bits). In most cases that
> we encounter, the desire for 387 FP is to gain the extra precision.
> Since very few CPUs support quad precision natively, the most
> obvious solution is a double-double library.

I also thought double-double would be sufficient, until last year, when
I ran into having to code an application for which it was unsuitable.
The problem with double-double was not in the number of bits in the
significand (106) -- the problem was that the exponent range of double
precision was not large enough. I believe the double-double type even
reduces the available exponent range from that of a double type. In any
case I ended up using the MPFR library for this application. IEEE quad
precision would be more suitable, since it expands the exponent range.

--
Krishna

The Cephes math library supports extended and 128b floats
https://netlib.org/cephes/128bdoc.html

But gnu libquadmath had been good enough for me until now.

On 18/04/2024 7:11 pm, albert@spenarnc.xs4all.nl wrote:
> ...
> Interesting read. The reverse polish calculator from HP was a
> resounding success, with great profits.
> Kahan contributed to this.
> It was killed by the bean counters at HP.
> There was huge demand but they refused to expand production.
> Then the calculator died because it simply was not available.

We had quite a few HP instruments but then the bottom fell out of
the test instrument market (and maintenance generally). It was
embarrassing to see the once great HP pursue the consumer PC market.

Stephen Pelc <stephen@vfxforth.com> writes:
>The manual (gasp) documents how to change the default FP package.
>
>Changing the default pack also changes the system call interfaces to
>match.

It is great that now an FP package is available by default, rather
than having to load some obscurely-named file from an
installation-specific path. E.g. appbench-1.4 contains a file
setup/vfx.fth which contains the lines:

1 cells 4 = [if]
include /usr/share/doc/VfxForth/Lib/Ndp387.fth
[then]

The result is that when I just tested various systems on the
appbench-1.4 suite, vfx64 worked (for four of the benchmarks), and I
won't spoil my positive message by reporting on what did not work.

- anton
--
M. Anton Ertl http://www.complang.tuwien.ac.at/anton/home.html
comp.lang.forth FAQs: http://www.complang.tuwien.ac.at/forth/faq/toc.html
New standard: https://forth-standard.org/
EuroForth 2023: https://euro.theforth.net/2023

mhx@iae.nl (mhx) writes:
>minforth wrote:
>> But I think the main advantage lies in the possibility of parallel and/or
>> vectorized execution.
>
>I have not yet seen algorithms where that would bring something.

Matrix multiplication is an easy case. I have also done a version of
Jon Bentley's greedy TSP program that benefitted from SSE and AVX; I
had to use assembly language to do this, however; see the thread
starting at <2016Nov14.164726@mips.complang.tuwien.ac.at>.

OTOH, yesterday I saw what gcc did for the inner loop of the bubble
benchmark from the Stanford integer benchmarks:

while ( i<top ) {

if ( sortlist[i] > sortlist[i+1] ) {
j = sortlist[i];
sortlist[i] = sortlist[i+1];
sortlist[i+1] = j;
};
i=i+1;
};

top=top-1;
};

gcc-12.2 -O1 produces straighforward scalar code, gcc-12.2 -O3 wants
to use SIMD instructions:

gcc -01 gcc -O3
1c: add $0x4,%rax c0: movq (%rax),%xmm0
cmp %rsi,%rax add $0x1,%edx
je 35 pshufd $0xe5,%xmm0,%xmm1
25: mov (%rax),%edx movd %xmm0,%edi
mov 0x4(%rax),%ecx movd %xmm1,%ecx
cmp %ecx,%edx cmp %ecx,%edi
jle 1c jle e1
mov %ecx,(%rax) pshufd $0xe1,%xmm0,%xmm0
mov %edx,0x4(%rax) movq %xmm0,(%rax)
jmp 1c e1: add $0x4,%rax
35: cmp %r8d,%edx
jl c0

The version produced by gcc -O3 is almost three times slower on a
Skylake than the one by gcc -O1 and is actually slower than several
Forth systems, including gforth-fast. I think that the reason is that
the movq towards the end stores two items, and the movq at the start
of the next iteration loads one of these item, i.e., there is partial
overlap between the store and the load. In this case the hardware
takes a slow path, which means that the slowdown is much bigger than
the instruction count suggests.

- anton
--
M. Anton Ertl http://www.complang.tuwien.ac.at/anton/home.html
comp.lang.forth FAQs: http://www.complang.tuwien.ac.at/forth/faq/toc.html
New standard: https://forth-standard.org/
EuroForth 2023: https://euro.theforth.net/2023

albert@spenarnc.xs4all.nl writes:
>I cut the same corners with ciforth. However I think this
>cannot be compliant with the IEEE requirement of the standard?

Which IEEE requirement of which standard?

- anton
--
M. Anton Ertl http://www.complang.tuwien.ac.at/anton/home.html
comp.lang.forth FAQs: http://www.complang.tuwien.ac.at/forth/faq/toc.html
New standard: https://forth-standard.org/
EuroForth 2023: https://euro.theforth.net/2023

anton@mips.complang.tuwien.ac.at (Anton Ertl) writes:
> The version produced by gcc -O3 is almost three times slower on a
> Skylake than the one by gcc -O1 and is actually slower than several
> Forth systems, including gforth-fast.

Wow, that seems worth a bug report. How about gcc -O2 ?

Paul Rubin <no.email@nospam.invalid> writes:
>anton@mips.complang.tuwien.ac.at (Anton Ertl) writes:
>> The version produced by gcc -O3 is almost three times slower on a
>> Skylake than the one by gcc -O1 and is actually slower than several
>> Forth systems, including gforth-fast.
>
>Wow, that seems worth a bug report.

At least the compiled code is working as the programmer intended. gcc
people regularly resolve bug reports as "INVALID" where earlier gcc
versions work as intended and later gcc versions do not. And then
there are the cases where they just do nothing, as on PR93811.

>How about gcc -O2 ?

I have not tried it.

If you want to measure it and/or report this as a bug, you can find
the source code in <http://www.complang.tuwien.ac.at/forth/bench.zip>
in the directory c-manual.

- anton
--
M. Anton Ertl http://www.complang.tuwien.ac.at/anton/home.html
comp.lang.forth FAQs: http://www.complang.tuwien.ac.at/forth/faq/toc.html
New standard: https://forth-standard.org/
EuroForth 2023: https://euro.theforth.net/2023

anton@mips.complang.tuwien.ac.at (Anton Ertl) writes:
>OTOH, yesterday I saw what gcc did for the inner loop of the bubble
>benchmark from the Stanford integer benchmarks:
>
> while ( i<top ) {
>
> if ( sortlist[i] > sortlist[i+1] ) {
> j = sortlist[i];
> sortlist[i] = sortlist[i+1];
> sortlist[i+1] = j;
> };
> i=i+1;
> };
>
> top=top-1;
> };
>
>gcc-12.2 -O1 produces straighforward scalar code, gcc-12.2 -O3 wants
>to use SIMD instructions:
>
> gcc -01 gcc -O3
>1c: add $0x4,%rax c0: movq (%rax),%xmm0
> cmp %rsi,%rax add $0x1,%edx
> je 35 pshufd $0xe5,%xmm0,%xmm1
>25: mov (%rax),%edx movd %xmm0,%edi
> mov 0x4(%rax),%ecx movd %xmm1,%ecx
> cmp %ecx,%edx cmp %ecx,%edi
> jle 1c jle e1
> mov %ecx,(%rax) pshufd $0xe1,%xmm0,%xmm0
> mov %edx,0x4(%rax) movq %xmm0,(%rax)
> jmp 1c e1: add $0x4,%rax
>35: cmp %r8d,%edx
> jl c0
>
>The version produced by gcc -O3 is almost three times slower on a
>Skylake than the one by gcc -O1 and is actually slower than several
>Forth systems, including gforth-fast. I think that the reason is that
>the movq towards the end stores two items, and the movq at the start
>of the next iteration loads one of these item, i.e., there is partial
>overlap between the store and the load. In this case the hardware
>takes a slow path, which means that the slowdown is much bigger than
>the instruction count suggests.

I was curious if a more recent Intel core had improved on that (and
maybe such a more recent Intel core was targeted by the "optimization"
that caused the slowdown), so I measured it on a P-core of a Core
i3-1315U. The results are as follows:

O1/bubble O3/bubble
424,248,952 2,061,809,866 cpu_core/cycles/
1,536,825,253 1,986,035,580 cpu_core/instructions/

So, more than a factor of 4 on this microarchitecture.

The differences in the topdown analysis are also interesting:

O1
1,177,188,340 cpu_core/topdown-retiring/ # 46.1% Retiring
279,332,826 cpu_core/topdown-bad-spec/ # 10.9% Bad Speculation
778,141,445 cpu_core/topdown-fe-bound/ # 30.5% Frontend Bound
319,237,516 cpu_core/topdown-be-bound/ # 12.5% Backend Bound
0 cpu_core/topdown-heavy-ops/ # 0.0% Heavy Operations
269,356,654 cpu_core/topdown-br-mispredict/ # 10.5% Branch Mispredict
269,356,654 cpu_core/topdown-fetch-lat/ # 10.5% Fetch Latency
59,857,034 cpu_core/topdown-mem-bound/ # 2.3% Memory Bound

O3
1,599,831,263 cpu_core/topdown-retiring/ # 12.9% Retiring
630,236,558 cpu_core/topdown-bad-spec/ # 5.1% Bad Speculation
533,277,087 cpu_core/topdown-fe-bound/ # 4.3% Frontend Bound
9,598,987,583 cpu_core/topdown-be-bound/ # 77.6% Backend Bound
280,169 cpu_core/topdown-heavy-ops/ # 0.0% Heavy Operations
630,236,558 cpu_core/topdown-br-mispredict/ # 5.1% Branch Mispredict
193,918,941 cpu_core/topdown-fetch-lat/ # 1.6% Fetch Latency
5,623,649,291 cpu_core/topdown-mem-bound/ # 45.5% Memory Bound

- anton
--
M. Anton Ertl http://www.complang.tuwien.ac.at/anton/home.html
comp.lang.forth FAQs: http://www.complang.tuwien.ac.at/forth/faq/toc.html
New standard: https://forth-standard.org/
EuroForth 2023: https://euro.theforth.net/2023

In article <2024Apr20.182321@mips.complang.tuwien.ac.at>,
Anton Ertl <anton@mips.complang.tuwien.ac.at> wrote:
>albert@spenarnc.xs4all.nl writes:
>>I cut the same corners with ciforth. However I think this
>>cannot be compliant with the IEEE requirement of the standard?
>
>Which IEEE requirement of which standard?

The ISO 9x talks about

64-bit IEEE double-precision number
32-bit IEEE double-precision number

IEEE floating point number as defined in ANSI/IEEE Standard
754-1985

80 bit '87 is not compliant with either of these, I think.

>- anton
--
Don't praise the day before the evening. One swallow doesn't make spring.
You must not say "hey" before you have crossed the bridge. Don't sell the
hide of the bear until you shot it. Better one bird in the hand than ten in
the air. First gain is a cat purring. - the Wise from Antrim -

albert@spenarnc.xs4all.nl writes:
>In article <2024Apr20.182321@mips.complang.tuwien.ac.at>,
>Anton Ertl <anton@mips.complang.tuwien.ac.at> wrote:
>>albert@spenarnc.xs4all.nl writes:
>>>I cut the same corners with ciforth. However I think this
>>>cannot be compliant with the IEEE requirement of the standard?
>>
>>Which IEEE requirement of which standard?
>
>The ISO 9x talks about

Whatever "ISO 9x" may be: Why do you worry about it and mention it
here?

> 64-bit IEEE double-precision number
> 32-bit IEEE double-precision number
>
>IEEE floating point number as defined in ANSI/IEEE Standard
>754-1985

So?

>80 bit '87 is not compliant with either of these, I think.

So?

- anton
--
M. Anton Ertl http://www.complang.tuwien.ac.at/anton/home.html
comp.lang.forth FAQs: http://www.complang.tuwien.ac.at/forth/faq/toc.html
New standard: https://forth-standard.org/
EuroForth 2023: https://euro.theforth.net/2023

On 21/04/2024 1:23 am, Anton Ertl wrote:
> Stephen Pelc <stephen@vfxforth.com> writes:
>> The manual (gasp) documents how to change the default FP package.
>>
>> Changing the default pack also changes the system call interfaces to
>> match.
>
> It is great that now an FP package is available by default, rather
> than having to load some obscurely-named file from an
> installation-specific path. E.g. appbench-1.4 contains a file
> setup/vfx.fth which contains the lines:
>
> 1 cells 4 = [if]
> include /usr/share/doc/VfxForth/Lib/Ndp387.fth
> [then]
>
> The result is that when I just tested various systems on the
> appbench-1.4 suite, vfx64 worked (for four of the benchmarks), and I
> won't spoil my positive message by reporting on what did not work.

What didn't work?

When loading an alternate f/p package I found it removes parts of the
compiler unrelated to fp. Hopefully this can be rectified. IIRC this
problem didn't exist with the previous setup where fp had to be explicitly
loaded.

On 22 Apr 2024 at 13:05:10 CEST, "dxf" <dxforth@gmail.com> wrote:

> When loading an alternate f/p package I found it removes parts of the
> compiler unrelated to fp. Hopefully this can be rectified. IIRC this
> problem didn't exist with the previous setup where fp had to be explicitly
> loaded.

Please send me more details by email.

Stephen
--
Stephen Pelc, stephen@vfxforth.com
MicroProcessor Engineering, Ltd. - More Real, Less Time
133 Hill Lane, Southampton SO15 5AF, England
tel: +44 (0)78 0390 3612, +34 649 662 974
http://www.mpeforth.com
MPE website
http://www.vfxforth.com/downloads/VfxCommunity/
downloads