Groups | Search | Server Info | Keyboard shortcuts | Login | Register [http] [https] [nntp] [nntps]

Groups > comp.lang.forth > #59451 > unrolled thread

matrix multiplication

Back to article view | Back to comp.lang.forth

Contents

  matrix multiplication krishna.myneni@ccreweb.org - 2017-05-17 15:13 -0700
    Re: matrix multiplication stephen@mpeforth.com (Stephen Pelc) - 2017-05-17 22:31 +0000
      Re: matrix multiplication krishna.myneni@ccreweb.org - 2017-05-17 15:39 -0700
        Re: matrix multiplication rickman <gnuarm@gmail.com> - 2017-05-18 03:22 -0400
          Re: matrix multiplication krishna.myneni@ccreweb.org - 2017-05-18 00:29 -0700
            Re: matrix multiplication rickman <gnuarm@gmail.com> - 2017-05-18 03:48 -0400
            Re: matrix multiplication anton@mips.complang.tuwien.ac.at (Anton Ertl) - 2017-05-18 08:42 +0000
    Re: matrix multiplication krishna.myneni@ccreweb.org - 2017-05-18 00:25 -0700
      Re: matrix multiplication anton@mips.complang.tuwien.ac.at (Anton Ertl) - 2017-05-19 15:13 +0000
        Re: matrix multiplication anton@mips.complang.tuwien.ac.at (Anton Ertl) - 2017-05-19 17:26 +0000
          Re: matrix multiplication Alex <alex@rivadpm.com> - 2017-05-19 18:56 +0100
          Re: matrix multiplication Paul Rubin <no.email@nospam.invalid> - 2017-05-19 21:48 -0700
          Re: matrix multiplication anton@mips.complang.tuwien.ac.at (Anton Ertl) - 2017-05-20 10:21 +0000
            Re: matrix multiplication krishna.myneni@ccreweb.org - 2017-05-20 06:12 -0700
            Re: matrix multiplication krishna.myneni@ccreweb.org - 2017-05-20 09:47 -0700
              Re: matrix multiplication Paul Rubin <no.email@nospam.invalid> - 2017-05-20 15:31 -0700
                Re: matrix multiplication krishna.myneni@ccreweb.org - 2017-05-21 04:46 -0700
                  Re: matrix multiplication mhx@iae.nl - 2017-05-21 13:41 -0700
                    Re: matrix multiplication krishna.myneni@ccreweb.org - 2017-05-21 15:32 -0700
              Re: matrix multiplication anton@mips.complang.tuwien.ac.at (Anton Ertl) - 2017-05-21 07:19 +0000
                Re: matrix multiplication krishna.myneni@ccreweb.org - 2017-05-21 05:10 -0700
                  Re: matrix multiplication anton@mips.complang.tuwien.ac.at (Anton Ertl) - 2017-05-21 14:29 +0000
                    Re: matrix multiplication krishna.myneni@ccreweb.org - 2017-05-21 09:17 -0700
        Re: matrix multiplication stephen@mpeforth.com (Stephen Pelc) - 2017-05-19 20:54 +0000
        Re: matrix multiplication krishna.myneni@ccreweb.org - 2017-05-19 14:33 -0700
          Re: matrix multiplication krishna.myneni@ccreweb.org - 2017-05-19 16:12 -0700
          Re: matrix multiplication anton@mips.complang.tuwien.ac.at (Anton Ertl) - 2017-05-20 14:48 +0000
      Re: matrix multiplication krishna.myneni@ccreweb.org - 2017-05-19 18:38 -0700

Page 1 of 2  [1] 2  Next page →

#59451 — matrix multiplication

Since there is an interesting discussion going on about matrix multiplication on integrated stack systems (allowed by Forth-94, but not by Forth-2012) vs separate fp stack systems, I spent an hour writing a full matrix multiplier for an integrated stack system, source below. I don't dispute the efficiency gained by a separate fp stack system, but I'm not terribly concerned by it either. The code below adheres to Forth-94, with the exception of the word A@ , which is simply @ in a standard system. With the factoring I used, it is clear that the inner-loop word DF_MUL_R1C2 is the bottleneck, and should be coded in assembler. It will be a useful exercise to do this and compare the execution times for all-Forth source vs mixed Forth and Forth-assembler words (the latter will vary from system to system). Note that there's no fancy algorithm here -- it is straightforward matrix multiplication.


Krishna

\ mmul.4th
\
\ Double precision floating point matrix multiplication
\ for an integrated data/fp stack.
\
\ Krishna Myneni, Creative Consulting for Research
\ and Education
\
\ Usage: a1 a2 a3 nr1 nc1 nc2 df_mmul
\
\ Requires:
\   ans-words.4th
\
\ Revisions:
\   2017-05-17  km; created.
\
\ Notes:
\   0. Matrix data is assumed to be stored in row order
\
\   1. Only the word DF_MMUL is specific for an integrated
\      data/fp stack. Other words work with separate fp
\      stack.

variable nc1
variable nc2
variable a1
variable a2

\ Convert row of a1 and column of a2 to
\ corresponding addresses
: df_r1c2>a1a2 ( row1 col2 -- arow1 acol2 )
    dfloats a2 a@ + >r
    nc1 @ dfloats * a1 a@ + r>   \ -- arow1 acol2
;

\ Multiply row of a1 with col of a2, element by element,
\ and accumulate the sum.
: df_mul_r1c2 ( row1 col2 -- rsum )
    df_r1c2>a1a2
    2>r 0e 2r>  \ rsum a1 a2
    nc1 @ 0 DO  
      2dup 2>r >r f@ r> f@ f* f+
      2r> nc2 @ dfloats + >r
      dfloat+ r>
    LOOP
    2drop
;

\ Multiply two double-precision matrices with data beginning at
\ a1 and a2, and store at a3. Proper memory allocation is
\ assumed, as are the dimensions for a2, i.e. nr2 = nc1 is
\ assumed. This word assumes an integrated data/fp stack.
: df_mmul ( a1 a2 a3 nr1 nc1 nc2 -- )
    nc2 ! nc1 ! 2>r a2 ! a1 !
    2r>  \ a3 nr1
    0 DO
      nc2 @ 0 DO
        J I df_mul_r1c2 2 pick f!
        dfloat+
      LOOP
    LOOP
    drop
;


\ Test code
\
\ Requires the additional files:
\   ttester.4th
\   fsl/fsl-util.4th

0 [IF]
base @
decimal

\ Allot and initialize three 2x2 matrices

2 2 double matrix a{{
2 2 double matrix b{{
2 2 double matrix c{{

t{ 1e a{{ 0 0 }} f! ->  }t
t{ 2e a{{ 0 1 }} f! ->  }t
t{ 3e a{{ 1 0 }} f! ->  }t
t{ 4e a{{ 1 1 }} f! ->  }t

t{ 5e b{{ 0 0 }} f! ->  }t
t{ 6e b{{ 0 1 }} f! ->  }t
t{ 7e b{{ 1 0 }} f! ->  }t
t{ 8e b{{ 1 1 }} f! ->  }t

TESTING df_mmul
t{ a{{ 0 0 }} b{{ 0 0 }} c{{ 0 0 }} 2 2 2 df_mmul -> }t
t{ c{{ 0 0 }} f@  ->  19e }t
t{ c{{ 0 1 }} f@  ->  22e }t
t{ c{{ 1 0 }} f@  ->  43e }t
t{ c{{ 1 1 }} f@  ->  50e }t

base !
[THEN]

[toc] | [next] | [standalone]

#59452

On Wed, 17 May 2017 15:13:49 -0700 (PDT), krishna.myneni@ccreweb.org
wrote:

>Since there is an interesting discussion going on about matrix multiplicati=
>on on integrated stack systems (allowed by Forth-94, but not by Forth-2012)=

As we keep trying to tell people, if there is a combined stack the
system has an environmental dependency. This means that you have to
document the dependency.

The effect of the wording changes is that a separate stack is the
default case.

Stephen

-- 
Stephen Pelc, stephenXXX@mpeforth.com
MicroProcessor Engineering Ltd - More Real, Less Time
133 Hill Lane, Southampton SO15 5AF, England
tel: +44 (0)23 8063 1441, fax: +44 (0)23 8033 9691
web: http://www.mpeforth.com - free VFX Forth downloads

[toc] | [prev] | [next] | [standalone]

#59453

On Wednesday, May 17, 2017 at 5:30:53 PM UTC-5, Stephen Pelc wrote:
> On Wed, 17 May 2017 15:13:49 -0700 (PDT), krishna....org
> wrote:
> 
> >Since there is an interesting discussion going on about matrix multiplicati=
> >on on integrated stack systems (allowed by Forth-94, but not by Forth-2012)=
> 
> As we keep trying to tell people, if there is a combined stack the
> system has an environmental dependency. This means that you have to
> document the dependency.

Wording is always difficult to get right. In your above sentence, I'm sure you mean, "the application has an environmental dependency," not "the system". Thanks for the clarification, though. I'll stop saying an integrated stack is "not allowed" by Forth-2012.

> 
> The effect of the wording changes is that a separate stack is the
> default case.
> 
> Stephen
> 

Krishna

[toc] | [prev] | [next] | [standalone]

#59464

On 5/17/2017 6:39 PM, krishna.myneni@ccreweb.org wrote:
> On Wednesday, May 17, 2017 at 5:30:53 PM UTC-5, Stephen Pelc wrote:
>> On Wed, 17 May 2017 15:13:49 -0700 (PDT), krishna....org
>> wrote:
>>
>>> Since there is an interesting discussion going on about matrix multiplicati=
>>> on on integrated stack systems (allowed by Forth-94, but not by Forth-2012)=
>>
>> As we keep trying to tell people, if there is a combined stack the
>> system has an environmental dependency. This means that you have to
>> document the dependency.
>
> Wording is always difficult to get right. In your above sentence, I'm sure you mean, "the application has an environmental dependency," not "the system". Thanks for the clarification, though. I'll stop saying an integrated stack is "not allowed" by Forth-2012.

It wouldn't be the program.  Your program is written to run on either 
type of system.  He should have said the system has an "environmental 
restriction".

-- 

Rick C

[toc] | [prev] | [next] | [standalone]

#59467

On Thursday, May 18, 2017 at 2:22:39 AM UTC-5, rickman wrote:
...
> It wouldn't be the program.  Your program is written to run on either 
> type of system.  ...

By program or application, I mean the Forth source, not the Forth environment. In the case of mmul.4th, it is written for an integrated data/fp stack. The word df_mmul will not work properly on a Forth system with a separate fp stack. It's a trivial mod to alter it to work with an fp_stack.

Krishna

[toc] | [prev] | [next] | [standalone]

#59469

On 5/18/2017 3:29 AM, krishna.myneni@ccreweb.org wrote:
> On Thursday, May 18, 2017 at 2:22:39 AM UTC-5, rickman wrote:
> ...
>> It wouldn't be the program.  Your program is written to run on either
>> type of system.  ...
>
> By program or application, I mean the Forth source, not the Forth environment. In the case of mmul.4th, it is written for an integrated data/fp stack. The word df_mmul will not work properly on a Forth system with a separate fp stack. It's a trivial mod to alter it to work with an fp_stack.

I gave some consideration once to the issues in writing a program to run 
on a system with an integrated stack.  Because of the potential 
difference in size of floating point data from the data stack width, I 
realized the only way to do so would be to work with floating point and 
other data entirely separately removing one from the top of stack before 
working with the other.  Trying to use stack manipulation instructions 
with both types of data is not possible if you don't know how many data 
stack cells the floating point data will use.

If the two types of data are essentially not mixed on the stack, the 
code should run on both types of systems.

How did you structure your stack usage that precludes working on a 
system with a separate stack?  I assume you must have made assumptions 
about the number of cells floating point data will use on the data stack?

-- 

Rick C

[toc] | [prev] | [next] | [standalone]

#59474

krishna.myneni@ccreweb.org writes:
>In the case of mmul.4th, it is written for an integrated data/fp stack.=
> The word df_mmul will not work properly on a Forth system with a separate =
>fp stack.

Assuming this is the only environmental dependency, the proper label
is to label it as a standard program with an environmental dependency
on a combined data/FP stack.  Note that you have to declare this
dependency for both Forth-94 and Forth-2012, because only programs
that work also on separate-stack systems are compliant without
qualification with Forth-94.

Conversely, a program that requires a separate FP stack would have to
be labeled as having an environmental dependency in Forth-94, and
needs no such qualification in Forth-2012.

For systems, in Forth-94 a system can choose either a separate stack
or any kind of combined stack without needing any qualification.  In
Forth-2012, a separate-stack system needs no qualification, while a
combined-stack system has an "environmental restriction of keeping
floating-point numbers on the data stack".

However, it is extremely unlikely that you have a dependency only on a
combined stack; such a program normally also has a dependency on a
specific FP item size.  And indeed, the code contains

2 pick f!

There is an FP item on the top-of-data stack at the start of this
word, so which item the PICK picks depends on the number of cells of
the FP item.

- 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: http://www.forth200x.org/forth200x.html
   EuroForth 2017: http://www.euroforth.org/ef17/

[toc] | [prev] | [next] | [standalone]

#59465

On Wednesday, May 17, 2017 at 5:13:50 PM UTC-5, krishna...@ccreweb.org wrote:
> ... I spent an hour writing a full matrix multiplier for an integrated stack system, source below. ...

The row offsets may be precomputed for the two matrices to be multiplied. The revised Forth source is given below. For multiplying two 100 x 100 double-precision floating point matrices (and storing the results in a third), I obtain an execution time of 96 ms, under kforth-fast, for the pure Forth code on a 1 GHz AMD Athlon-64 system (four-core system, but using a single process):

$ kforth-fast ans-words ttester fsl/fsl-util fsl/horner fsl/extras/noise mmul

FSL-UTIL          V1.19f         29 April   2010   EFC F2DUP is redefined

HORNER            V1.5d          19 February  2012   EFC 
NOISE             V1.0e         12 March     2013   KM TESTING df_mmul
 ok

100 100 double matrix a{{
 ok
100 100 double matrix b{{
 ok
100 100 double matrix c{{
 ok
100 0 do 100 0 do ran0 a{{ j i }} f! loop loop
 ok
100 0 do 100 0 do ran0 b{{ j i }} f! loop loop
 ok
20 0 do us2@ a{{ 0 0 }} b{{ 0 0 }} c{{ 0 0 }} 100 100 100 df_mmul us2@ 2swap d- loop
 ok

 ok
20 0 do d. cr loop
96074 
100124 
95517 
95947 
96245 
96718 
96511 
97235 
95988 
95881 
95531 
95901 
96243 
110905 
108721 
116355 
117974 
131328 
143829 
234237 
 ok

\ Elapsed times above are in microseconds.

-----------------------------------------------
\ mmul.4th
\
\ Double precision floating point matrix multiplication
\ for an integrated data/fp stack.
\
\ Krishna Myneni, Creative Consulting for Research
\ and Education
\
\ Usage: a1 a2 a3 nr1 nc1 nc2 df_mmul
\
\ Requires:
\   ans-words.4th
\
\ Revisions:
\   2017-05-17  km; created.
\   2017-05-18  km; precompute row offsets for a1 and a2.
\
\ Notes:
\   0. Matrix data is assumed to be stored in row order
\
\   1. Only the word DF_MMUL is specific for an integrated
\      data/fp stack. Other words work with separate fp
\      stack.

variable nc1
variable nc2
variable a1
variable a2
variable roffs1
variable roffs2

\ Convert row of a1 and column of a2 to
\ corresponding addresses
: df_r1c2>a1a2 ( row1 col2 -- arow1 acol2 )
    dfloats a2 a@ + >r
    roffs1 @ * a1 a@ + r>   \ -- arow1 acol2
;

\ Multiply row of a1 with col of a2, element by element,
\ and accumulate the sum.
: df_mul_r1c2 ( row1 col2 -- rsum )
    df_r1c2>a1a2
    2>r 0e 2r>  \ rsum a1 a2
    nc1 @ 0 DO  
      2dup 2>r >r f@ r> f@ f* f+
      2r> roffs2 @ + >r
      dfloat+ r>
    LOOP
    2drop
;

\ Multiply two double-precision matrices with data beginning at
\ a1 and a2, and store at a3. Proper memory allocation is
\ assumed, as are the dimensions for a2, i.e. nr2 = nc1 is
\ assumed. This word assumes an integrated data/fp stack.
: df_mmul ( a1 a2 a3 nr1 nc1 nc2 -- )
    nc2 ! nc1 ! 2>r a2 ! a1 !

    \ offsets to next row for a1 and a2
    nc1 @ dfloats roffs1 !
    nc2 @ dfloats roffs2 !

    2r>  \ a3 nr1
    0 DO
      nc2 @ 0 DO
        J I df_mul_r1c2 2 pick f!
        dfloat+
      LOOP
    LOOP
    drop
;


\ Test code
\
\ Requires the additional files:
\   ttester.4th
\   fsl/fsl-util.4th

true [IF]
base @
decimal

\ Allot and initialize three 2x2 matrices

2 2 double matrix a{{
2 2 double matrix b{{
2 2 double matrix c{{

t{ 1e a{{ 0 0 }} f! ->  }t
t{ 2e a{{ 0 1 }} f! ->  }t
t{ 3e a{{ 1 0 }} f! ->  }t
t{ 4e a{{ 1 1 }} f! ->  }t

t{ 5e b{{ 0 0 }} f! ->  }t
t{ 6e b{{ 0 1 }} f! ->  }t
t{ 7e b{{ 1 0 }} f! ->  }t
t{ 8e b{{ 1 1 }} f! ->  }t

TESTING df_mmul
t{ a{{ 0 0 }} b{{ 0 0 }} c{{ 0 0 }} 2 2 2 df_mmul -> }t
t{ c{{ 0 0 }} f@  ->  19e }t
t{ c{{ 0 1 }} f@  ->  22e }t
t{ c{{ 1 0 }} f@  ->  43e }t
t{ c{{ 1 1 }} f@  ->  50e }t

base !
[THEN]

[toc] | [prev] | [next] | [standalone]

#59559

krishna.myneni@ccreweb.org writes:
>On Wednesday, May 17, 2017 at 5:13:50 PM UTC-5, krishna...@ccreweb.org wrot=
>e:
>> ... I spent an hour writing a full matrix multiplier for an integrated st=
>ack system, source below. ...
>
>The row offsets may be precomputed for the two matrices to be multiplied. T=
>he revised Forth source is given below.

I took this code, replaced occurences of A@ with @ and replaced

>        J I df_mul_r1c2 2 pick f!
>        dfloat+

with

        J I rot >r df_mul_r1c2 r@ f!
        r> dfloat+

and the result seems to be working; I did not find this replacement
trivial, BTW.  I also produced a version of the inner loop that works
only on separate-stack systems (i.e., Forth-2012):

Forth-94                            Forth-2012
nc1 @ 0 DO                          nc1 @ 0 DO        
  2dup 2>r >r f@ r> f@ f* f+          2dup f@ f@ f* f+
  2r> roffs2 @ + >r                   roffs2 @ + >r   
  dfloat+ r>                          dfloat+ r>      
LOOP                                LOOP              

Results for a 500x500 matrix multiply on a Core i5-6600K:

Forth-94 Forth-2012 Forth-2012 Forth-94 Forth-94 
Myneni  Myneni/Ertl     Ertl   Rees/Ertl   pahihu 
109966M     107831M    1046M       1202M    1517M vfxlin Ndp387.fth
  2986M       2260M    2220M       2945M    2764M gforth-fast
  4028M       4076M    4844M       5067M    5017M SwiftForth HW stack, no WAIT

I don't know what in this code slows down VFX so much.  It does not
appear to be the classical code-close-to-written-data problem (the
inner loop does not write data close to the code, and the number of
icache misses (2.8M) is not high enough to explain such a big number
of cycles).  The code generated by VFX (see below) for the inner loop
looks inconspicuous to me.

SwiftForth benefits from having only 5 (instead of 8 for the other
variants) FP words (i.e., call-return pairs) in the inner loop.

Here's the code that VFX generates for the inner loop

( 080C6E60    DB2B )                  FLD       TBYTE 0 [EBX]
( 080C6E62    8B5500 )                MOV       EDX, [EBP]
( 080C6E65    DB2A )                  FLD       TBYTE 0 [EDX]
( 080C6E67    DEC9 )                  FMULP     ST(1), ST
( 080C6E69    DEC1 )                  FADDP     ST(1), ST
( 080C6E6B    031D70340A08 )          ADD       EBX, [080A3470]
( 080C6E71    53 )                    PUSH      EBX
( 080C6E72    8B5D00 )                MOV       EBX, [EBP]
( 080C6E75    83C308 )                ADD       EBX, 08
( 080C6E78    5A )                    POP       EDX
( 080C6E79    83042401 )              ADD       [ESP], 01
( 080C6E7D    8344240401 )            ADD       [ESP+04], 01
( 080C6E82    895D00 )                MOV       [EBP], EBX
( 080C6E85    8BDA )                  MOV       EBX, EDX
( 080C6E87    71D7 )                  JNO       080C6E60

Finally, Gforth also contains a primitive V* that fits the inner loop
of your matrix multiplication:

v* ( f_addr1 nstride1 f_addr2 nstride2 ucount -- r )

I implemented a variant of your inner word that uses V*:

: df_mul_r1c2 ( row1 col2 -- rsum )
    df_r1c2>a1a2 [ 1 floats ] literal swap roffs2 @ nc1 @ v* ;

Let's compare the version that uses V* with the version that uses FAXPY:

cycles
  faxpy       v*
   316M     562M  500x500
19.168M 120.710M 2000x2000

For the 500x500 case, the V* version is slower because the latency of
the F+ is a lower bound for the performance per iteration of V*,
because each F+ depends on the result of the F+ of the previous
iteration; by contrast, in faxpy the floating-point operations of each
iteration are independent of other iterations.

For the 2000x2000 case, the V* version is slower because it has no
spatial locality when accessing the column, so for each element of the
column, we see a cache miss, while the FAXPY version accesses
everything in a spatially local way; hmm, swapping the middle and the
outer loop around might help, by providing at least temporal locality
for the column.

- 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: http://www.forth200x.org/forth200x.html
   EuroForth 2017: http://www.euroforth.org/ef17/

[toc] | [prev] | [next] | [standalone]

#59560

anton@mips.complang.tuwien.ac.at (Anton Ertl) writes:
>Results for a 500x500 matrix multiply on a Core i5-6600K:
>
>Forth-94 Forth-2012 Forth-2012 Forth-94 Forth-94 
>Myneni  Myneni/Ertl     Ertl   Rees/Ertl   pahihu 
>109966M     107831M    1046M       1202M    1517M vfxlin Ndp387.fth
>  2986M       2260M    2220M       2945M    2764M gforth-fast
>  4028M       4076M    4844M       5067M    5017M SwiftForth HW stack, no WAIT
>
>I don't know what in this code slows down VFX so much.  It does not
>appear to be the classical code-close-to-written-data problem (the
>inner loop does not write data close to the code, and the number of
>icache misses (2.8M) is not high enough to explain such a big number
>of cycles).  The code generated by VFX (see below) for the inner loop
>looks inconspicuous to me.

I looked at it with "perf record -e cycles" and "perf report", and got
the following for the inner loop:

       ( 080C6E60    FF7500 )       PUSH      [EBP]
 0.26% ( 080C6E63    53 )           PUSH      EBX
       ( 080C6E64    53 )           PUSH      EBX
       ( 080C6E65    8B5D00 )       MOV       EBX, [EBP]
       ( 080C6E68    DB2B )         FLD       TBYTE 0 [EBX]
 0.85% ( 080C6E6A    5B )           POP       EBX
       ( 080C6E6B    DB2B )         FLD       TBYTE 0 [EBX]
 1.75% ( 080C6E6D    DEC9 )         FMULP     ST(1), ST
52.46% ( 080C6E6F    DEC1 )         FADDP     ST(1), ST
42.25% ( 080C6E71    5B )           POP       EBX
 1.11% ( 080C6E72    5A )           POP       EDX
       ( 080C6E73    031D70340A08 ) ADD       EBX, [080A3470]
 0.23% ( 080C6E79    53 )           PUSH      EBX
 0.23% ( 080C6E7A    83C208 )       ADD       EDX, 08
       ( 080C6E7D    5B )           POP       EBX
 0.32% ( 080C6E7E    83042401 )     ADD       [ESP], 01
 0.12% ( 080C6E82    8344240401 )   ADD       [ESP+04], 01
 0.20% ( 080C6E87    895500 )       MOV       [EBP], EDX
       ( 080C6E8A    71D4 )         JNO       080C6E60

Everything else had 0.08% or fewer cycles.  The locations associated
with the events are not very precise, but in this case they point to
the FADDP.  It's still unclear to me what the problem is.  One thing
that one might suspect is denormal numbers, but they cannot occur with
this input (only integers converted to floats).

- 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: http://www.forth200x.org/forth200x.html
   EuroForth 2017: http://www.euroforth.org/ef17/

[toc] | [prev] | [next] | [standalone]

#59561

On 5/19/2017 18:26, Anton Ertl wrote:
> anton@mips.complang.tuwien.ac.at (Anton Ertl) writes:
>> Results for a 500x500 matrix multiply on a Core i5-6600K:
>>
>> Forth-94 Forth-2012 Forth-2012 Forth-94 Forth-94
>> Myneni  Myneni/Ertl     Ertl   Rees/Ertl   pahihu
>> 109966M     107831M    1046M       1202M    1517M vfxlin Ndp387.fth
>>   2986M       2260M    2220M       2945M    2764M gforth-fast
>>   4028M       4076M    4844M       5067M    5017M SwiftForth HW stack, no WAIT
>>
>> I don't know what in this code slows down VFX so much.  It does not
>> appear to be the classical code-close-to-written-data problem (the
>> inner loop does not write data close to the code, and the number of
>> icache misses (2.8M) is not high enough to explain such a big number
>> of cycles).  The code generated by VFX (see below) for the inner loop
>> looks inconspicuous to me.
> 
> I looked at it with "perf record -e cycles" and "perf report", and got
> the following for the inner loop:
> 
>         ( 080C6E60    FF7500 )       PUSH      [EBP]
>   0.26% ( 080C6E63    53 )           PUSH      EBX
>         ( 080C6E64    53 )           PUSH      EBX
>         ( 080C6E65    8B5D00 )       MOV       EBX, [EBP]
>         ( 080C6E68    DB2B )         FLD       TBYTE 0 [EBX]
>   0.85% ( 080C6E6A    5B )           POP       EBX
>         ( 080C6E6B    DB2B )         FLD       TBYTE 0 [EBX]
>   1.75% ( 080C6E6D    DEC9 )         FMULP     ST(1), ST
> 52.46% ( 080C6E6F    DEC1 )         FADDP     ST(1), ST
> 42.25% ( 080C6E71    5B )           POP       EBX
>   1.11% ( 080C6E72    5A )           POP       EDX
>         ( 080C6E73    031D70340A08 ) ADD       EBX, [080A3470]
>   0.23% ( 080C6E79    53 )           PUSH      EBX
>   0.23% ( 080C6E7A    83C208 )       ADD       EDX, 08
>         ( 080C6E7D    5B )           POP       EBX
>   0.32% ( 080C6E7E    83042401 )     ADD       [ESP], 01
>   0.12% ( 080C6E82    8344240401 )   ADD       [ESP+04], 01
>   0.20% ( 080C6E87    895500 )       MOV       [EBP], EDX
>         ( 080C6E8A    71D4 )         JNO       080C6E60
> 
> Everything else had 0.08% or fewer cycles.  The locations associated
> with the events are not very precise, but in this case they point to
> the FADDP.  It's still unclear to me what the problem is.  One thing
> that one might suspect is denormal numbers, but they cannot occur with
> this input (only integers converted to floats).
> 
> - anton
> 

Have you tried inserting a NOP between the FMULP and FADDP? The FADP 
falls over a 16 byte paragraph boundary.

-- 
Alex

[toc] | [prev] | [next] | [standalone]

#59574

anton@mips.complang.tuwien.ac.at (Anton Ertl) writes:
>>I don't know what in this code slows down VFX so much.

32 bit vs 64 bit mode?

[toc] | [prev] | [next] | [standalone]

#59586

anton@mips.complang.tuwien.ac.at (Anton Ertl) writes:
>anton@mips.complang.tuwien.ac.at (Anton Ertl) writes:
>>Results for a 500x500 matrix multiply on a Core i5-6600K:
>>
>>Forth-94 Forth-2012 Forth-2012 Forth-94 Forth-94 
>>Myneni  Myneni/Ertl     Ertl   Rees/Ertl   pahihu 
>>109966M     107831M    1046M       1202M    1517M vfxlin Ndp387.fth
>>  2986M       2260M    2220M       2945M    2764M gforth-fast
>>  4028M       4076M    4844M       5067M    5017M SwiftForth HW stack, no WAIT
>>
>>I don't know what in this code slows down VFX so much.

Found it.  The problem was that my benchmark uses FLOATS/FLOAT+, while
Krishna Myneni's code uses DFLOATS/DFLOAT+.  This worked ok on Gforth
and SwiftForth where floats and dfloats have the same size, but on VFX
the sizes are different, so the FP items processed by the matrix
multiplication code were not the intended FP items.  My guess is that
a significant number of them were denormal numbers and other
exceptional cases that are not processed fast by the FPU.

I fixed it by converting Krishna Myneni's code to use FLOATS and
FLOAT+, and voila:

Forth-94 Forth-2012 F-2012 Forth-2012 Forth-94 Forth-94 
Myneni  Myneni/Ertl Myneni+   Ertl   Rees/Ertl   pahihu 
  1770M       1252M  1170M   1046M       1202M    1517M vfxlin Ndp387.fth
  2986M       2260M  2152M   2220M       2945M    2764M gforth-fast
  4028M       4076M  4077M   4844M       5067M    5017M SwiftForth ...

The Myneni+ variant adds Stephen Pelc' suggestion of using SWAP FLOAT+
SWAP instead of >R FLOAT+ R>.

- 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: http://www.forth200x.org/forth200x.html
   EuroForth 2017: http://www.euroforth.org/ef17/

[toc] | [prev] | [next] | [standalone]

#59595

On Saturday, May 20, 2017 at 5:53:59 AM UTC-5, Anton Ertl wrote:
...
> 
> Forth-94 Forth-2012 F-2012 Forth-2012 Forth-94 Forth-94 
> Myneni  Myneni/Ertl Myneni+   Ertl   Rees/Ertl   pahihu 
>   1770M       1252M  1170M   1046M       1202M    1517M vfxlin Ndp387.fth
>   2986M       2260M  2152M   2220M       2945M    2764M gforth-fast
>   4028M       4076M  4077M   4844M       5067M    5017M SwiftForth ...
> 
> The Myneni+ variant adds Stephen Pelc' suggestion of using SWAP FLOAT+
> SWAP instead of >R FLOAT+ R>.
> 

Stephen's suggestion to replace >R DFLOAT+ R> with SWAP DFLOAT+ SWAP improves the performance in kforth-fast also. For multiplying 100x100 double-precision matrices, the execution time is reduced by 16%. Implementation of SWAP is critical. Last time I checked, Intel's XCHG instruction to exchange the contents of a register and memory is slower than using MOV.

Krishna

[toc] | [prev] | [next] | [standalone]

#59601

On Saturday, May 20, 2017 at 5:53:59 AM UTC-5, Anton Ertl wrote:
...
> Forth-94 Forth-2012 F-2012 Forth-2012 Forth-94 Forth-94 
> Myneni  Myneni/Ertl Myneni+   Ertl   Rees/Ertl   pahihu 
>   1770M       1252M  1170M   1046M       1202M    1517M vfxlin Ndp387.fth
>   2986M       2260M  2152M   2220M       2945M    2764M gforth-fast
>   4028M       4076M  4077M   4844M       5067M    5017M SwiftForth ...
> 
> The Myneni+ variant adds Stephen Pelc' suggestion of using SWAP FLOAT+
> SWAP instead of >R FLOAT+ R>.
> 

A next step could be to do the computations in parallel on a multi-core system, following the example I posted earlier:

ftp://ccreweb.org/software/kforth/examples/parallel-mm.4th

Do gforth, vfxlin, and swift forth provide fork and mmap syscalls, or equivalents? We should be able to realize a speedup of ~n times by distributing the computing over n cores. Matrix multiplication is an easy computation to parallelize.

I'll try writing a parallel calculating version of the mmul.4th module for testing. Eventually, the inner loop word, df_mul_r1c2 will be rewritten in assembler.

Krishna

[toc] | [prev] | [next] | [standalone]

#59614

krishna.myneni@ccreweb.org writes:
> A next step could be to do the computations in parallel on a
> multi-core system...  Do gforth, vfxlin, and swift forth provide fork
> and mmap syscalls, or equivalents?

Before that I'd suggest generating some vector code to use the Intel
multimedia instructions or on-chip graphics processor, and after that
use a GPU.

[toc] | [prev] | [next] | [standalone]

#59641

On Saturday, May 20, 2017 at 5:31:06 PM UTC-5, Paul Rubin wrote:
> krishna...org writes:
> > A next step could be to do the computations in parallel on a
> > multi-core system...  Do gforth, vfxlin, and swift forth provide fork
> > and mmap syscalls, or equivalents?
> 
> Before that I'd suggest generating some vector code to use the Intel
> multimedia instructions or on-chip graphics processor, and after that
> use a GPU.


Well, I'm just going after the low-hanging fruit, which doesn't require considerable effort. I already know how to use a forked process to obtain ~2x speed increase for a matrix multiplication computation, so it wouldn't be much work for me to do that. However, use of SSE2 (or later) instructions would certainly improve the efficiency of each process (by up to 4x?). I just haven't played with this yet. As for GPUs, I don't think they're relevant yet for double-precision fp calculations.

Krishna

[toc] | [prev] | [next] | [standalone]

#59697

On Sunday, May 21, 2017 at 1:46:57 PM UTC+2, krishna...@ccreweb.org wrote:
> On Saturday, May 20, 2017 at 5:31:06 PM UTC-5, Paul Rubin wrote:
> > krishna...org writes:
> > > A next step could be to do the computations in parallel on a
> > > multi-core system...  Do gforth, vfxlin, and swift forth provide fork
> > > and mmap syscalls, or equivalents?
> > 
> > Before that I'd suggest generating some vector code to use the Intel
> > multimedia instructions or on-chip graphics processor, and after that
> > use a GPU.
[..] 
 
> As for GPUs, I don't think they're relevant yet for
> double-precision fp calculations.
 
Nvidia GeForce GTX 970: 109 double-precision Gflops for ~200USD$
Radeon R9 390X: 790G dpflops for ~400USD$
The single-precision performance is 10 times higher.

-marcel

[toc] | [prev] | [next] | [standalone]

#59702

On Sunday, May 21, 2017 at 3:42:00 PM UTC-5, m...@iae.nl wrote:
> On Sunday, May 21, 2017 at 1:46:57 PM UTC+2, krishna...org wrote:
> > On Saturday, May 20, 2017 at 5:31:06 PM UTC-5, Paul Rubin wrote:
> > > krishna...org writes:
> > > > A next step could be to do the computations in parallel on a
> > > > multi-core system...  Do gforth, vfxlin, and swift forth provide fork
> > > > and mmap syscalls, or equivalents?
> > > 
> > > Before that I'd suggest generating some vector code to use the Intel
> > > multimedia instructions or on-chip graphics processor, and after that
> > > use a GPU.
> [..] 
>  
> > As for GPUs, I don't think they're relevant yet for
> > double-precision fp calculations.
>  
> Nvidia GeForce GTX 970: 109 double-precision Gflops for ~200USD$
> Radeon R9 390X: 790G dpflops for ~400USD$
> The single-precision performance is 10 times higher.

Looks like my information is outdated. Thanks.

Krishna

[toc] | [prev] | [next] | [standalone]

#59637

krishna.myneni@ccreweb.org writes:
>A next step could be to do the computations in parallel on a multi-core system, following the example I posted earlier:
>
>ftp://ccreweb.org/software/kforth/examples/parallel-mm.4th
>
>Do gforth, vfxlin, and swift forth provide fork and mmap syscalls, or equivalents?

There is multi-threading support in these systems (for Gforth in the
development version) that you could use; then you also don't need to
share the memory with mmap.

However, my next step would be to implement cache blocking.  While the
slowdown from going beyond the cache is (to me) surprisingly low for
the FAXPY-based version, that is on a system that is otherwise idle
(interestingly both on my PC and on the Odroid C2).  If the system
makes more use of its multiple cores, and especially if you want to
have several instances of matrix multiply running on them,using less
main memory bandwidth by using cache blocking should help.

Another thing that I would (in general) do before using multiple cores
is vectorization.  Unfortunately, we don't have vector support in
Forth, not even in a system-specific way, yet (AFAIK).

And if I am really interested in multiplying a lot of matrices
quickly, I would use a library specialized for that, i.e., ATLAS, or
one provided by the CPU manufacturer.

> We should be able to realize a speedup of ~n times by distributing the computing over n cores. Matrix multiplication is an easy computation to parallelize.

But you realize this speedup only if your computation does not become
memory-bandwidth 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: http://www.forth200x.org/forth200x.html
   EuroForth 2017: http://www.euroforth.org/ef17/

[toc] | [prev] | [next] | [standalone]

Page 1 of 2  [1] 2  Next page →

Back to top | Article view | comp.lang.forth

csiph-web