Fix BUG-013: generalize create_chunks off concrete Ferrite.Grid

.github/workflows/CI.yml

@@ -47,6 +47,8 @@ jobs:
     name: Documentation
     runs-on: ubuntu-latest
     timeout-minutes: 30
+    env:
+      JULIA_NUM_THREADS: 2
     steps:
         - uses: actions/checkout@v7
         - uses: julia-actions/setup-julia@v3

docs/src/literate_tutorials/iga.jl

@@ -81,6 +81,23 @@ apply!(a, ch)
 work!(assembler, buffer; a=a)
 apply!(r, lh, 0.0);
 
+# ## Threaded IGA assembly
+# `setup_domainbuffer` also supports `threading=true` for `BezierGrid`, since
+# `FerriteAssembly` no longer requires a concrete `Ferrite.Grid` for chunk creation.
+# We request 2 tasks explicitly (independent of `Threads.nthreads()`) so the
+# chunk-splitting is exercised even in a single-threaded Julia session, and
+# compare against the sequential result above.
+threaded_buffer = setup_domainbuffer(domain; threading=true, num_tasks=2)
+K_threaded = allocate_matrix(dh)
+r_threaded = zeros(ndofs(dh))
+assembler_threaded = start_assemble(K_threaded, r_threaded)
+work!(assembler_threaded, threaded_buffer; a=a)
+apply!(r_threaded, lh, 0.0);
+
+using Test                #hide
+@test K_threaded ≈ K      #hide
+@test r_threaded ≈ r      #hide
+
 # before solving it,
 apply_zero!(K, r, ch)
 a .-= K\r
@@ -111,9 +128,8 @@ FerriteIGA.VTKIGAFile("plate_with_hole.vtu", grid) do vtk
     write_solution(vtk, dh, a)
 end
 
-using Test                                      #src
 # @test sum(norm, σ_nodes) ≈ 3087.2447327126742 #src
-@test norm(norm.(qe.data)) ≈ 679.3207411544098  #src
+@test norm(norm.(qe.data)) ≈ 679.3207411544098  #hide
 
 #md # ## [Plain program](@id iga_plain_program)
 #md #

src/Multithreading/TaskChunks.jl

@@ -81,7 +81,7 @@ function split_in_chunks(set::Vector{T}; num_tasks = Threads.nthreads()) where T
 end
 
 # If chunks already given, check that they match the intersected set
-function create_chunks(::Grid, intersected_set::Vector{I}, chunks::Vector{Vector{Vector{I}}}) where I
+function create_chunks(::Ferrite.AbstractGrid, intersected_set::Vector{I}, chunks::Vector{Vector{Vector{I}}}) where I
     chunk_set = sizehint!(Set{I}(), length(intersected_set))
     full_set = Set(intersected_set)
     for chunk_vector in chunks
@@ -100,22 +100,22 @@ function create_chunks(::Grid, intersected_set::Vector{I}, chunks::Vector{Vector
     return chunks
 end
 
-# Colors given with same type as set 
-function create_chunks(g::Grid, intersected_set::Vector{I}, colors::Vector{Vector{I}}) where I
+# Colors given with same type as set
+function create_chunks(g::Ferrite.AbstractGrid, intersected_set::Vector{I}, colors::Vector{Vector{I}}) where I
     colors_intersect = map(sort! ∘ collect ∘ Base.Fix1(intersect, intersected_set), colors)
     chunks = [split_in_chunks(set) for set in colors_intersect]
     return create_chunks(g, intersected_set, chunks)
 end
 # Colors given with different type as set
-function create_chunks(g::Grid, intersected_set::Vector{FacetIndex}, colors::Vector{Vector{Int}})
+function create_chunks(g::Ferrite.AbstractGrid, intersected_set::Vector{FacetIndex}, colors::Vector{Vector{Int}})
     cellset = first.(intersected_set)
     colors_intersect = map(sort! ∘ collect ∘ Base.Fix1(intersect, cellset), colors)
     chunks = [convert_chunk(split_in_chunks(set), intersected_set) for set in colors_intersect]
     return create_chunks(g, intersected_set, chunks) # Pass through to check correctness
 end
 
 # Colors not given
-function create_chunks(grid::Grid, intersected_set::Vector, ::Nothing)
+function create_chunks(grid::Ferrite.AbstractGrid, intersected_set::Vector, ::Nothing)
     makecellset(v::Vector{Int}) = Set(v)
     makecellset(v::Vector) = Set(first.(v))
     return create_chunks(grid, intersected_set, create_coloring(grid, makecellset(intersected_set)))

test/threading_utils.jl

@@ -59,6 +59,31 @@ end
     end
 end
 
+# `create_chunks` previously dispatched on concrete `Ferrite.Grid`, so any other
+# `Ferrite.AbstractGrid` implementation (e.g. FerriteIGA.jl's `BezierGrid`) could not
+# be used with `threading=true`, even though `create_coloring` and the generic
+# `AbstractGrid` accessors it relies on only require `.cells` and `.nodes` fields.
+struct DummyGrid{C,N} <: Ferrite.AbstractGrid{2}
+    cells::Vector{C}
+    nodes::Vector{N}
+end
+
+@testset "create_chunks with non-Grid AbstractGrid (BUG-013)" begin
+    grid = generate_grid(Quadrilateral, (4, 4))
+    dummygrid = DummyGrid(grid.cells, grid.nodes)
+    cellset = collect(1:getncells(grid))
+
+    # Automatic coloring path (colors_or_chunks = nothing)
+    chunks_grid = FerriteAssembly.create_chunks(grid, cellset, nothing)
+    chunks_dummy = FerriteAssembly.create_chunks(dummygrid, cellset, nothing)
+    @test Set(Iterators.flatten(Iterators.flatten(chunks_dummy))) == Set(cellset)
+    @test [sort!(collect(Iterators.flatten(c))) for c in chunks_dummy] == [sort!(collect(Iterators.flatten(c))) for c in chunks_grid]
+
+    # User-supplied chunks path (doesn't use the grid argument beyond dispatch)
+    chunks = [[cellset[1:8], cellset[9:16]], [cellset[17:end]]]
+    @test FerriteAssembly.create_chunks(dummygrid, cellset, chunks) == chunks
+end
+
 @testset "work! with empty custom chunks (PR89)" begin
     # End-to-end reproduction: a domain whose user-supplied chunks contain
     # empty sub-chunks must still visit every cell during threaded `work!`.