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 # SVP64 Power ISA version
  10
  11 Summary: 9 instructions (see below, 8 instructions) 5 of which are 64-bit
  12 for a total of 14 "words".
  13
  14 ```
  15     # r5: n count; r6: x ptr; r7: y ptr; fp1: a
  16     1  mtctr 5                # move n to CTR
  17     2  addi r10,r6,0          # copy y-ptr into r10 (y')
  18     3  .L2
  19     4  setvl MAXVL=32,VL=CTR  # actually VL=MIN(MAXVL,CTR)
  20     5  sv.lfdup *32,8(6)      # load x into fp32-63, inc x
  21     6  sv.lfdup *64,8(7)      # load y into fp64-95, inc y
  22     7  sv.fmadd *64,*64,1,*32 # (*y) = (*y) * (*x) + fp1
  23     8  sv.stfdup *64,8(10)    # store at y-copy, inc y'
  24     9  sv.bc/ctr .L2          # decr CTR by VL, jump !zero
  25     10 blr                    # return
  26 ```
  27
  28 A refinement, reducing 1 instruction and register port usage.
  29 Relies on post-increment, relies on no overlap between x and y
  30 in memory, and critically relies on y overwrite.
  31
  32 ```
  33     # r5: n count; r6: x ptr; r7: y ptr; fp1: a
  34     1  mtctr 5                # move n to CTR
  35     2  .L2
  36     3  setvl MAXVL=32,VL=CTR  # actually VL=MIN(MAXVL,CTR)
  37     4  sv.lfdup *32,8(6)      # load x into fp32-63, incr x
  38     5  sv.lfd *64,8(7)        # load y into fp64-95, NO INC
  39     6  sv.fmadd *64,*64,1,*32 # (*y) = (*y) * (*x) + fp1
  40     7  sv.stfdup *64,8(7)     # store at y, incr y
  41     8  sv.bc/ctr .L2          # decr CTR by VL, jump !zero
  42     9  blr                    # return
  43 ```
  44
  45 # RVV version
  46
  47 Summary: 12 instructions, 7 32-bit and 5 16-bit for a total of 9.5 "words"
  48
  49 ```
  50     # a0 is n, a1 is pointer to x[0], a2 is pointer to y[0], fa0 is a
  51       li t0, 2<<25
  52       vsetdcfg t0             # enable 2 64b Fl.Pt. registers
  53     loop:
  54       setvl  t0, a0           # vl = t0 = min(mvl, n)
  55       vld    v0, a1           # load vector x
  56       c.slli   t1, t0, 3      # t1 = vl * 8 (in bytes)
  57       vld    v1, a2           # load vector y
  58       c.add    a1, a1, t1     # increment pointer to x by vl*8
  59       vfmadd v1, v0, fa0, v1  # v1 += v0 * fa0 (y = a * x + y)
  60       c.sub    a0, a0, t0     # n -= vl (t0)
  61       vst    v1, a2           # store Y
  62       c.add    a2, a2, t1     # increment pointer to y by vl*8
  63       c.bnez   a0, loop       # repeat if n != 0
  64       c.ret                   # return
  65 ```
  66
  67 # SVE Version
  68
  69 Summary: 12 instructions, all 32-bit, for a total of 12 "words"
  70
  71 ```
  72     1  // x0 = &x[0], x1 = &y[0], x2 = &a, x3 = &n
  73     2  daxpy_:
  74     3  ldrswx3, [x3]                     // x3=*n
  75     4  movx4, #0                         // x4=i=0
  76     5  whilelt p0.d, x4, x3              // p0=while(i++<n)
  77     6  ld1rdz0.d, p0/z, [x2]             // p0:z0=bcast(*a)
  78     7  .loop:
  79     8  ld1d z1.d, p0/z, [x0, x4, lsl #3] // p0:z1=x[i]
  80     9  ld1d z2.d, p0/z, [x1, x4, lsl #3] // p0:z2=y[i]
  81     10 fmla z2.d, p0/m, z1.d, z0.d       // p0?z2+=x[i]*a
  82     11 st1d z2.d, p0, [x1, x4, lsl #3]   // p0?y[i]=z2
  83     12 incd x4                           // i+=(VL/64)
  84     13 .latch:
  85     14 whilelt p0.d, x4, x3             // p0=while(i++<n)
  86     15 b.first .loop                    // more to do?
  87     16 ret
  88 ```