From: lkcl <lkcl@web>
Date: Sat, 15 Apr 2023 12:42:38 +0000 (+0100)
Subject: (no commit message)
X-Git-Tag: opf_rfc_ls009_v1~79
X-Git-Url: https://git.libre-soc.org/?a=commitdiff_plain;h=b1598707cb45898286d3eb0a05d86e95a2a40e97;p=libreriscv.git

---

diff --git a/openpower/sv/rfc/ls009.mdwn b/openpower/sv/rfc/ls009.mdwn
index c41016ace..b25d94ce7 100644
--- a/openpower/sv/rfc/ls009.mdwn
+++ b/openpower/sv/rfc/ls009.mdwn
@@ -1758,5 +1758,176 @@ if __name__ == '__main__':
     demo()
 ```
 
+## REMAP FFT pseudocode
+
+```
+# a "yield" version of the REMAP algorithm, for FFT Tukey-Cooley schedules
+# original code for the FFT Tukey-Cooley schedule:
+# https://www.nayuki.io/res/free-small-fft-in-multiple-languages/fft.py
+"""
+    # Radix-2 decimation-in-time FFT (real, not complex)
+    size = 2
+    while size <= n:
+        halfsize = size // 2
+        tablestep = n // size
+        for i in range(0, n, size):
+            k = 0
+            for j in range(i, i + halfsize):
+                jh = j+halfsize
+                jl = j
+                temp1 = vec[jh] * exptable[k]
+                temp2 = vec[jl]
+                vec[jh] = temp2 - temp1
+                vec[jl] = temp2 + temp1
+                k += tablestep
+        size *= 2
+"""
+
+# python "yield" can be iterated. use this to make it clear how
+# the indices are generated by using natural-looking nested loops
+def iterate_butterfly_indices(SVSHAPE):
+    # get indices to iterate over, in the required order
+    n = SVSHAPE.lims[0]
+    stride = SVSHAPE.lims[2] # stride-multiplier on reg access
+    # creating lists of indices to iterate over in each dimension
+    # has to be done dynamically, because it depends on the size
+    # first, the size-based loop (which can be done statically)
+    x_r = []
+    size = 2
+    while size <= n:
+        x_r.append(size)
+        size *= 2
+    # invert order if requested
+    if SVSHAPE.invxyz[0]: x_r.reverse()
+
+    if len(x_r) == 0:
+        return
+
+    # start an infinite (wrapping) loop
+    skip = 0
+    while True:
+        for size in x_r:           # loop over 3rd order dimension (size)
+            x_end = size == x_r[-1]
+            # y_r schedule depends on size
+            halfsize = size // 2
+            tablestep = n // size
+            y_r = []
+            for i in range(0, n, size):
+                y_r.append(i)
+            # invert if requested
+            if SVSHAPE.invxyz[1]: y_r.reverse()
+            for i in y_r:       # loop over 2nd order dimension
+                y_end = i == y_r[-1]
+                k_r = []
+                j_r = []
+                k = 0
+                for j in range(i, i+halfsize):
+                    k_r.append(k)
+                    j_r.append(j)
+                    k += tablestep
+                # invert if requested
+                if SVSHAPE.invxyz[2]: k_r.reverse()
+                if SVSHAPE.invxyz[2]: j_r.reverse()
+                for j, k in zip(j_r, k_r):   # loop over 1st order dimension
+                    z_end = j == j_r[-1]
+                    # now depending on MODE return the index
+                    if SVSHAPE.skip == 0b00:
+                        result = j              # for vec[j]
+                    elif SVSHAPE.skip == 0b01:
+                        result = j + halfsize   # for vec[j+halfsize]
+                    elif SVSHAPE.skip == 0b10:
+                        result = k              # for exptable[k]
+
+                    loopends = (z_end |
+                               ((y_end and z_end)<<1) |
+                                ((y_end and x_end and z_end)<<2))
+
+                    yield (result * stride) + SVSHAPE.offset, loopends
+
+def demo():
+    # set the dimension sizes here
+    xdim = 8
+    ydim = 0 # not needed
+    zdim = 1 # stride must be set to 1
+
+    # set total. err don't know how to calculate how many there are...
+    # do it manually for now
+    VL = 0
+    size = 2
+    n = xdim
+    while size <= n:
+        halfsize = size // 2
+        tablestep = n // size
+        for i in range(0, n, size):
+            for j in range(i, i + halfsize):
+                VL += 1
+        size *= 2
+
+    # set up an SVSHAPE
+    class SVSHAPE:
+        pass
+    # j schedule
+    SVSHAPE0 = SVSHAPE()
+    SVSHAPE0.lims = [xdim, ydim, zdim]
+    SVSHAPE0.order = [0,1,2]  # experiment with different permutations, here
+    SVSHAPE0.mode = 0b01
+    SVSHAPE0.skip = 0b00
+    SVSHAPE0.offset = 0       # experiment with different offset, here
+    SVSHAPE0.invxyz = [0,0,0] # inversion if desired
+    # j+halfstep schedule
+    SVSHAPE1 = SVSHAPE()
+    SVSHAPE1.lims = [xdim, ydim, zdim]
+    SVSHAPE1.order = [0,1,2]  # experiment with different permutations, here
+    SVSHAPE0.mode = 0b01
+    SVSHAPE1.skip = 0b01
+    SVSHAPE1.offset = 0       # experiment with different offset, here
+    SVSHAPE1.invxyz = [0,0,0] # inversion if desired
+    # k schedule
+    SVSHAPE2 = SVSHAPE()
+    SVSHAPE2.lims = [xdim, ydim, zdim]
+    SVSHAPE2.order = [0,1,2]  # experiment with different permutations, here
+    SVSHAPE0.mode = 0b01
+    SVSHAPE2.skip = 0b10
+    SVSHAPE2.offset = 0       # experiment with different offset, here
+    SVSHAPE2.invxyz = [0,0,0] # inversion if desired
+
+    # enumerate over the iterator function, getting new indices
+    schedule = []
+    for idx, (jl, jh, k) in enumerate(zip(iterate_butterfly_indices(SVSHAPE0),
+                                          iterate_butterfly_indices(SVSHAPE1),
+                                          iterate_butterfly_indices(SVSHAPE2))):
+        if idx >= VL:
+            break
+        schedule.append((jl, jh, k))
+
+    # ok now pretty-print the results, with some debug output
+    size = 2
+    idx = 0
+    while size <= n:
+        halfsize = size // 2
+        tablestep = n // size
+        print ("size %d halfsize %d tablestep %d" % \
+                (size, halfsize, tablestep))
+        for i in range(0, n, size):
+            prefix = "i %d\t" % i
+            k = 0
+            for j in range(i, i + halfsize):
+                (jl, je), (jh, he), (ks, ke) = schedule[idx]
+                print ("  %-3d\t%s j=%-2d jh=%-2d k=%-2d -> "
+                        "j[jl=%-2d] j[jh=%-2d] ex[k=%d]" % \
+                                (idx, prefix, j, j+halfsize, k,
+                                      jl, jh, ks,
+                                ),
+                                "end", bin(je)[2:], bin(je)[2:], bin(ke)[2:])
+                k += tablestep
+                idx += 1
+        size *= 2
+
+# run the demo
+if __name__ == '__main__':
+    demo()
+```
+
+
 
 [[!tag opf_rfc]]