simple_v_extension/daxpy_example.mdwn

   1 # c code
   2
   3 ```
   4     void daxpy(size_t n, double a, const double x[], double y[]) {
   5         for (size_t i = 0; i < n; i++)
   6             y[i] = a*x[i] + y[i];
   7     }
   8 ```
   9
  10 Summary
  11
  12 | ISA | total | loop | words | notes |
  13 |-----|-------|------|-------|-------|
  14 | SVP64 | 9 | 7 | 14 | 5 64-bit, 4 32-bit |
  15 | RVV | 13 | 11 | 9.5 | 7 32-bit, 5 16-bit |
  16 | SVE | 12 | 7 | 12 | all 32-bit |
  17
  18 # SVP64 Power ISA version
  19
  20 Relies on post-increment, relies on no overlap between x and y
  21 in memory, and critically relies on y overwrite. x is post-incremented
  22 when read, but y is post-incremented on write. Element-Strided
  23 ensures the Immediate (8) results in a contiguous LD (or store)
  24 despite RA being marked Scalar (*without* modifying RA, on `sv.lfd/els`).
  25 For `sv.lfdup`, RA is Scalar so that only one
  26 LD/ST Update "wins": the last write to RA is the address for
  27 the next block.
  28
  29 ```
  30     # r5: n count; r6: x ptr; r7: y ptr; fp1: a
  31     1  addi r3,r7,0           # return result
  32     2  mtctr 5                # move n to CTR
  33     3  .L2
  34     4  setvl MAXVL=32,VL=CTR  # actually VL=MIN(MAXVL,CTR)
  35     5  sv.lfdup/els *32,8(6)  # load x into fp32-63, incr x
  36     6  sv.lfd/els *64,8(7)    # load y into fp64-95, NO INC
  37     7  sv.fmadd *64,*64,1,*32 # (*y) = (*y) * (*x) + a
  38     8  sv.stfdup/els *64,8(7) # store at y, incr y
  39     9  sv.bc/ctr .L2          # decr CTR by VL, jump !zero
  40     10 blr                    # return
  41 ```
  42
  43 # RVV version
  44
  45
  46 ```
  47     # a0 is n, a1 is pointer to x[0], a2 is pointer to y[0], fa0 is a
  48       li t0, 2<<25
  49       vsetdcfg t0             # enable 2 64b Fl.Pt. registers
  50     loop:
  51       setvl  t0, a0           # vl = t0 = min(mvl, n)
  52       vld    v0, a1           # load vector x
  53       c.slli   t1, t0, 3      # t1 = vl * 8 (in bytes)
  54       vld    v1, a2           # load vector y
  55       c.add    a1, a1, t1     # increment pointer to x by vl*8
  56       vfmadd v1, v0, fa0, v1  # v1 += v0 * fa0 (y = a * x + y)
  57       c.sub    a0, a0, t0     # n -= vl (t0)
  58       vst    v1, a2           # store Y
  59       c.add    a2, a2, t1     # increment pointer to y by vl*8
  60       c.bnez   a0, loop       # repeat if n != 0
  61       c.ret                   # return
  62 ```
  63
  64 # SVE Version
  65
  66 ```
  67     1  // x0 = &x[0], x1 = &y[0], x2 = &a, x3 = &n
  68     2  daxpy_:
  69     3  ldrswx3, [x3]                     // x3=*n
  70     4  movx4, #0                         // x4=i=0
  71     5  whilelt p0.d, x4, x3              // p0=while(i++<n)
  72     6  ld1rdz0.d, p0/z, [x2]             // p0:z0=bcast(*a)
  73     7  .loop:
  74     8  ld1d z1.d, p0/z, [x0, x4, lsl #3] // p0:z1=x[i]
  75     9  ld1d z2.d, p0/z, [x1, x4, lsl #3] // p0:z2=y[i]
  76     10 fmla z2.d, p0/m, z1.d, z0.d       // p0?z2+=x[i]*a
  77     11 st1d z2.d, p0, [x1, x4, lsl #3]   // p0?y[i]=z2
  78     12 incd x4                           // i+=(VL/64)
  79     13 .latch:
  80     14 whilelt p0.d, x4, x3             // p0=while(i++<n)
  81     15 b.first .loop                    // more to do?
  82     16 ret
  83 ```