More on Lisp structures performance

This is a follow up on my earlier email about the performance difference between regular Lisp structures and custom packed structures serialised into a byte array. There were certain issues in how the benchmarking was performed (in short, the tests were too short and the CPU may have not gotten a chance to rev the clock speed up) and I was getting conflicting results. I have addressed that and did some more investigation on how to get the most out of the LW compiler when it comes to numerical throughput and there are some conclusions that I came to:

1. LispWorks structures can actually be quite fast in practice, good enough to make some computations acceptable for production code.
2. Custom serialised structs implemented on top of typed-aref constructs can be made even faster, if one is willing to put time into custom struct development.
3. … but even the latter still won’t match performance of native code generated by C (or SBCL for that matter) and that does not have anything to do with some magical obscure compiler quirks and tricks.

Let’s unpack it a bit.

To make long story short, here is the vanilla struct code that we were using:

(defstruct particle
  (id          0 :type (unsigned-byte 32))
  (last-col-id 0 :type (unsigned-byte 32))
  (x     0.0     :type single-float)
  (y     0.0     :type single-float)
  (vel-x 0.0     :type single-float)
  (vel-y 0.0     :type single-float)
  (mass  1.0     :type single-float)
  (size  1.0     :type single-float))

Above is the Lisp structure itself that we all know and (this is optional) love. This is a structure representing a 2D particle struct that contains its x and y coordinates, as well as its velocity along either axis, using which we can update the location of the particle. The rest of the fields are not used and are really just padding to make the task for the memory bus a bit more challenging and the structure look more real-life like.

However we need to keep in mind that here each field is represented by an immediate 64 bit word, essentially giving us a packed linear representation of the data. This is what allows us to retrieve the values quickly — we just offset into a memory location within the struct, rather than chase a pointer pointing to some remote location in the RAM. This is what we have to do for double floats, since while technically doubles do fit into 64 bits just fine (that is their native size, why wouldn’t they?), but lisp requires each value it owns to have a passport to identify itself, which in LW is an additional tag of 3 bits, which makes now a double not to fit into the field. But for other values it works just fine, check this diagram of an immediate single float layout within a struct field:


64-BIT LISPWORKS WORD (IMMEDIATE SINGLE FLOAT)
      =========================================

Bit:  63                               32 31            8 7      0
      +--------------------------------+----------------+--------+
      |                                |                |        |
      |   IEEE 754 Single Payload      |   Zero / Pad   |   TAG  |
      |                                |                |        |
      +--------------------------------+----------------+--------+
      |S| Exp (8) |   Mantissa (23)    | 0000 .... 0000 |11110111|
      +--------------------------------+----------------+--------+
       ^                                                 ^
       |                                                 |
  Float Value                                        Tag = 247
 (Upper 32 bits)                                    (Hex #xF7)

In LW an immediate float is represented by a 64 bit value that has its lower bits 111 (value of 7, possibly to identify an other-pointer value), plus a few more bits in the free area in the lower 32-bit half of the word that give us more specific information on what type of other-pointer this value is. The upper half contains the raw float that could be directly loaded into a single float register (those are commonly identified by Sn symbols, as in s1, s2, s3, etc).

We can confirm that this layout is correct by:

(format nil “#x~x" (sys:object-pointer 1.0))
“#x3F800000000000F7"

Have a look at the last byte, which is #xF7 confirming the suggested layout. Such a layout could be of interest for some further data packing hacks, for example, we could fit a 24 bit integer into the bits 32-8 to identify the square in the world grid where the particle is located (or use two 12 bit integers for that), or we could use that space as 3 x 8 bit bytes to store, for example, the colour of the particle. Whether LW would allow us to do it is a different question that I have not thought of yet, but it may be possible in theory.

Let us move to the next part of the code.

(declaim (inline particle-x particle-y particle-vel-x particle-vel-y 
                 (setf particle-x) (setf particle-y)))

Here we simply tell the compiler to write the assembly code necessary to access the fields of the structure into the body of the function that accesses it during the compilation stage, rather than make a full blown function call every time a value in the structutre is needed to be read, or written to.

(defun update-particles-fast (particles dt)
  ;; the usual declarations and optimisations
  (declare (type (simple-array particle (*)) particles)
           (type single-float dt)
           (optimize (speed 3) (safety 0) (debug 0) (float 0)))
  (dotimes (x 100)
    (dotimes (i +num-particles+)
      (let ((p (aref particles i)))
        (setf (particle-x p) (+ (particle-x p) 
                                (* (particle-vel-x p) dt)))
        (setf (particle-y p) (+ (particle-y p) 
                                (* (particle-vel-y p) dt)))))))

In this part of the code we are performing 100 batches of particle processing, each batch is 10M particles, making a test of computing and updating 1B structures. For reference, our galaxy contains about 100B stars, which means that a test that could aim to model the Milky Way would be practical to run on commodity hardware, although you would need a hefty amount of ram to fit all the stars in (and RAM is hard to come by these days due to the AI gobbling up all the world supply for this year).

Updating the stars' coordinates would be rather straight forward using our code above, however if we were to model their gravitational interactions as well, our test would fall apart quickly — we would need to calculate the interactions between each pair of the stars, which would be 100B^2 / 2. That is 0.5d22 in Common Lisp double float notation, a rather substantial figure, or 5,000,000,000,000,000,000,000 in human numbers in case it does not look impressive enough.

A bit too much of a test for a lappy on the couch while watching netflix, but we won’t worry about it now but will worry about it some other time. Onto the next part.

(defun test-structs nil
  (clean-down)                             ; warm cpu up
  (update-particles-fast *particles* +dt+) ; warm it even more!
  (format t "~2%2. Running FAST update (Speed 3, Safety 0, Inline):~%")
  (time (update-particles-fast *particles* +dt+))) ; warm enough now

The above is the code for the actual test, nothing to see here really, it is here for completeness. The only thing of interest is that in order to get representative figures from our benchmarking, we need to warm the CPU up first, so to speak, to shift it from the idle low clock speed mode into the full tilt high clock speed mode. From the testing it also turned out that running a GC pass by (clean-down) is normally sufficient for a warm-up but I put an additional particle update run without timing just to be completely clean.

So now that we have the book-keeping details out of the way, here are the testing results:

2. Running FAST update (Speed 3, Safety 0, Inline):
Timing the evaluation of (STRUCT-TEST::UPDATE-PARTICLES-FAST STRUCT-TEST::*PARTICLES* STRUCT-TEST::+DT+)

User time    =        1.974
System time  =        0.018
Elapsed time =        1.981
Allocation   = 162672 bytes
11 Page faults
GC time      =        0.000

We ran the 1B particle test in just a tad less than 2 seconds. That gives us a throughput of close to 500M particles per second, or 

(* 1.981 3.2d9 1d-9)

~6.34 cycles per struct on the CPU that runs at 3.2GHz. Not bad, given that a single struct pass requires 12 cycles to compute:

1. 4 memory reads.
2. 4 arithmetic operations.
3. 2 memory writes.
4. 2 cycles of loop overhead (increment and compare).

We essentially had a double throughput of the scalar performance we would expect from a CPU with the given clock speed. This is evidently where the parallel execution pipeline architecture in ARM processors lends us a helping hand, but the twist of the plot is — we have 4 of such threads on the processor, so why are seeing a double performance increase, rather than 4x? Seems like we should do better.

How much better we can do will be discussed in another email, but before that, does anyone believe we could break the 1B particles per second sonic boom barrier? While remaining in the scalar domain, of course (i.e. computing just one particle at a time). Take note of your guess and we’ll see if that was correct later.