Parallel step races: allocation counters, the stepping thread's index, stage clocks

aephysics/basics/module.ae

@@ -18,6 +18,7 @@
 // hold. The allocation counters the original keeps (b3Alloc, b3Free with
 // sizes) are kept so a leak is a number.
 import std.string
+import aephysics.native
 
 exports (
     BitSet, IdPool, HashSet, Ints, Buffer, Stack, Arena, BlockAllocator, BLOCK_SIZE,
@@ -57,22 +58,25 @@ var g_bytes_in_use = 0
 
 // Box3D's b3Alloc: zeroed, 16-byte aligned by the C runtime's own promise,
 // and counted, so a leak is a number a test reads.
+// The counters move atomically: the collide pass and finalize allocate from
+// several workers at once (a mesh contact's triangle cache, an arena's
+// overflow, the continuous pass's sensor hits).
 alloc(size: int) -> ptr {
     if size <= 0 { return null }
-    g_alloc_count = g_alloc_count + 1
-    g_bytes_in_use = g_bytes_in_use + size
+    native.atomic_add((&g_alloc_count) as ptr, 1)
+    native.atomic_add((&g_bytes_in_use) as ptr, size)
     return calloc(1, size)
 }
 
 free_bytes(block: ptr, size: int) {
     if block == null { return }
-    g_alloc_count = g_alloc_count - 1
-    g_bytes_in_use = g_bytes_in_use - size
+    native.atomic_add((&g_alloc_count) as ptr, 0 - 1)
+    native.atomic_add((&g_bytes_in_use) as ptr, 0 - size)
     free(block)
 }
 
-alloc_count() -> int { return g_alloc_count }
-bytes_in_use() -> int { return g_bytes_in_use }
+alloc_count() -> int { return native.atomic_load((&g_alloc_count) as ptr) }
+bytes_in_use() -> int { return native.atomic_load((&g_bytes_in_use) as ptr) }
 
 // --- bits (ctz.h) ----------------------------------------------------------
 

aephysics/parallel/module.ae

@@ -244,8 +244,15 @@ var g_for_handles: ptr = null      // *Task[MAX_WORKERS]
 parallel_for(s: *Scheduler, body: fn(int, int, int, ptr) -> void, item_count: int, min_range: int, context: ptr) -> int {
     if item_count <= 0 { return 0 }
     worker_count = s.worker_count
+    // The caller works as worker 0 -- the task contexts and every module's
+    // scratch are chosen by the thread's worker index -- and gets its own
+    // index back after: a host stepping a one-worker world from one of its
+    // pool's threads (index k) indexed the world's single task context at k.
+    caller_index = native.worker_index()
     if worker_count <= 1 {
+        native.set_worker_index(0)
         body(0, item_count, 0, context)
+        native.set_worker_index(caller_index)
         return 1
     }
     if g_for_tasks == null {
@@ -278,5 +285,6 @@ parallel_for(s: *Scheduler, body: fn(int, int, int, ptr) -> void, item_count: in
         finish(s, handles[i] as *Task)
         i = i + 1
     }
+    native.set_worker_index(caller_index)
     return block_count
 }

aephysics/physics_world/module.ae

@@ -133,6 +133,11 @@ is_valid(id: WorldId) -> bool { return dynamics.world_is_valid(id) }
 step(id: WorldId, time_step: float, sub_step_count: int) {
     w = world_of(id)
     if w == null { return }
+    // The step runs on the calling thread as worker 0 (its task context and
+    // scratch), whatever the thread's index in a host's pool; the index is
+    // given back at the end.
+    caller_index = native.worker_index()
+    native.set_worker_index(0)
     solver.stage_step()
     w.locked = true
     basics.buffer_clear(&w.body_move_events)
@@ -182,6 +187,7 @@ step(id: WorldId, time_step: float, sub_step_count: int) {
     dynamics.swap_end_events(w)
     w.profile.step = ms(os.now_monotonic_ns() - t0)
     w.locked = false
+    native.set_worker_index(caller_index)
 }
 
 ms(ns: long) -> float { return (ns as float) / 1000000.0 }

aephysics/solver/module.ae

@@ -62,30 +62,48 @@ stages_timed() -> bool {
     return g_stage_on == 1
 }
 
-stage_clock() -> float[] {
-    if g_stage_ms == null { g_stage_ms = basics.alloc(STAGE_COUNT * 8) }
-    return g_stage_ms as float[]
-}
-
+// A row of stage clocks per worker, so the workers timing their blocks at
+// once neither race to allocate it nor lose each other's sums; made on the
+// stepping thread (stage_step) and summed for the report.
 stage_now_ms() -> float { return (os.now_monotonic_ns() as float) / 1000000.0 }
 
 stage_note(kind: int, started: float) {
-    if kind < 0 || kind >= STAGE_COUNT { return }
-    clock = stage_clock()
+    if kind < 0 || kind >= STAGE_COUNT || g_stage_ms == null { return }
+    clock = native.scratch(g_stage_ms, STAGE_COUNT * 8) as float[]
     clock[kind] = clock[kind] + (stage_now_ms() - started)
 }
 
 // One more step over the same stages, for the averages.
 stage_step() {
-    if stages_timed() { g_stage_steps = g_stage_steps + 1 }
+    if stages_timed() {
+        if g_stage_ms == null { g_stage_ms = basics.alloc(native.MAX_WORKERS * STAGE_COUNT * 8) }
+        g_stage_steps = g_stage_steps + 1
+    }
+}
+
+// The workers' rows summed.
+stage_clock() -> float[] {
+    totals = basics.alloc(STAGE_COUNT * 8) as float[]
+    rows = g_stage_ms as float[]
+    worker = 0
+    while worker < native.MAX_WORKERS {
+        kind = 0
+        while kind < STAGE_COUNT {
+            totals[kind] = totals[kind] + rows[worker * STAGE_COUNT + kind]
+            kind = kind + 1
+        }
+        worker = worker + 1
+    }
+    return totals
 }
 
 // The averages, in the reference's own order and words.
 stage_report() {
-    if stages_timed() == false || g_stage_steps == 0 { return }
+    if stages_timed() == false || g_stage_steps == 0 || g_stage_ms == null { return }
     clock = stage_clock()
     n = g_stage_steps as float
     println("aephysics stages ms/step: prepare_joints ${clock[STAGE_PREPARE_JOINTS] / n} prepare_wide ${clock[STAGE_PREPARE_WIDE] / n} prepare_contacts ${clock[STAGE_PREPARE_CONTACTS] / n} integrate_v ${clock[STAGE_INTEGRATE_VELOCITIES] / n} warm ${clock[STAGE_WARM_START] / n} solve ${clock[STAGE_SOLVE] / n} integrate_p ${clock[STAGE_INTEGRATE_POSITIONS] / n} relax ${clock[STAGE_RELAX] / n} restitution ${clock[STAGE_RESTITUTION] / n} store ${clock[STAGE_STORE_WIDE] / n}")
+    basics.free_bytes(clock as ptr, STAGE_COUNT * 8)
 }
 
 const NULL_INDEX = 0 - 1

tests/test_physics_world.ae

@@ -20,6 +20,7 @@ import aephysics.distance
 import aephysics.broad_phase
 import aephysics.dynamics
 import aephysics.physics_world
+import aephysics.native
 import aephysics.joint_solver
 
 extern exit(code: int)
@@ -454,6 +455,58 @@ test_joint_springs_and_locks() {
     hull.destroy_hull(cube)
 }
 
+// The parallel step's bookkeeping: a four-worker world allocating from its
+// workers (every new mesh contact grows a triangle cache in the collide pass)
+// frees everything it counted; and a step from a thread that is not worker 0
+// of anything (a host pool's thread) runs as worker 0 and gives the index back.
+test_parallel_bookkeeping() {
+    before = basics.alloc_count()
+    wdef = dynamics.default_world_def()
+    wdef.worker_count = 4
+    world_id = dynamics.create_world(&wdef)
+    grid = mesh.create_grid_mesh(16, 16, 1.0, 1, true)
+    gdef = dynamics.default_body_def()
+    ground = dynamics.create_body(world_id, &gdef)
+    sdef = dynamics.default_shape_def()
+    dynamics.create_mesh_shape(ground, &sdef, grid, math.vec3_one())
+    bdef = dynamics.default_body_def()
+    bdef.kind = dynamics.BODY_DYNAMIC
+    i = 0
+    while i < 200 {
+        bdef.position = math.vec3(-7.5 + ((i % 15) as float), 0.6 + 0.01 * ((i / 15) as float), -7.0 + ((i / 15) as float))
+        ball = dynamics.create_body(world_id, &bdef)
+        dynamics.create_sphere_shape(ball, &sdef, manifold.sphere(math.vec3_zero(), 0.4))
+        i = i + 1
+    }
+    i = 0
+    while i < 20 {
+        physics_world.step(world_id, 1.0 / 60.0, 4)
+        i = i + 1
+    }
+    dynamics.destroy_world(world_id)
+    mesh.destroy_mesh(grid)
+    ensure("parallel: everything the workers counted freed (${basics.alloc_count() - before})", basics.alloc_count() == before)
+
+    wdef = dynamics.default_world_def()
+    world_id = dynamics.create_world(&wdef)
+    slab = hull.make_box_hull(4.0, 0.5, 4.0)
+    gdef.position = math.vec3(0.0, 0.0 - 0.5, 0.0)
+    ground = dynamics.create_body(world_id, &gdef)
+    dynamics.create_hull_shape(ground, &sdef, slab)
+    bdef.position = math.vec3(0.0, 0.45, 0.0)
+    ball = dynamics.create_body(world_id, &bdef)
+    dynamics.create_sphere_shape(ball, &sdef, manifold.sphere(math.vec3_zero(), 0.5))
+    native.set_worker_index(5)
+    physics_world.step(world_id, 1.0 / 60.0, 4)
+    physics_world.step(world_id, 1.0 / 60.0, 4)
+    after = native.worker_index()
+    native.set_worker_index(0)
+    ensure("parallel: a host thread's index given back (${after})", after == 5)
+    ensure("parallel: the ball rests on the slab", dynamics.world_contact_count(dynamics.get_world(world_id.index1 - 1)) == 1)
+    dynamics.destroy_world(world_id)
+    hull.destroy_hull(slab)
+}
+
 test_hit_events() {
     wdef = dynamics.default_world_def()
     wdef.hit_event_threshold = 1.0
@@ -973,6 +1026,7 @@ main() {
     test_rebuilt_static_pairs()
     test_restitution()
     test_joint_springs_and_locks()
+    test_parallel_bookkeeping()
     test_continuous_move_event()
     test_continuous_gravity()
     test_sensor_stepped()