Make close of ConstraintHandler 50x faster

[KnutAM]

While making #859 I discovered that close!(ch) was quite slow. Did a profile and turns out most time was spent on copy!(ch.free_dofs, setdiff(1:ndofs(ch.dh), ch.prescribed_dofs)). For the benchmarked test case this PR reduces the construction of a ConstraintHandler from about 2 s to 0.04 s.
(Take exact speedup with grain of salt since done on laptop using battery power)

using Ferrite
grid = generate_grid(Quadrilateral, (1000, 1000))

dh = DofHandler(grid)
ip = Lagrange{RefQuadrilateral,2}()
add!(dh, :u, ip^2)
add!(dh, :v, ip)
close!(dh)

function make_ch(dh, grid)
    ch = ConstraintHandler(dh)
    add!(ch, Dirichlet(:u, getfaceset(grid, "right"), Returns(zero(Vec{2}))))
    add!(ch, Dirichlet(:v, getfaceset(grid, "left"), Returns(0.0)))
    close!(ch)
    return ch
end

ch = make_ch(dh, grid)
@time make_ch(dh, grid)

# PR
# run 1: 0.037628 seconds (4.26 k allocations: 93.643 MiB, 27.37% gc time)
# run 2: 0.027672 seconds (4.26 k allocations: 93.643 MiB)
  
# master 
# run 1: 1.989513 seconds (4.28 k allocations: 528.475 MiB, 2.04% gc time)
# run 2: 1.974789 seconds (4.28 k allocations: 528.475 MiB)
  

[fredrikekre]

Nice, perhaps extract to a sorted_setdiff! helper?

[codecov-commenter]

Codecov Report

All modified and coverable lines are covered by tests ✅

Comparison is base (cfe3863) 93.28% compared to head (7961299) 93.29%.

❗ Your organization needs to install the Codecov GitHub app to enable full functionality.

Additional details and impacted files

@@            Coverage Diff             @@
##           master     #860      +/-   ##
==========================================
+ Coverage   93.28%   93.29%   +0.01%     
==========================================
  Files          36       36              
  Lines        5224     5235      +11     
==========================================
+ Hits         4873     4884      +11     
  Misses        351      351              

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

[KristofferC]

I wonder if there should exist a free_dofs variable computed at all? It seems wasteful to have a long array of numbers that is very close to the dense set of 1:ndofs. I wonder if is_free(x, dh) = !is_constrained(x, dh) is good enough for example.

[KnutAM]

If storing it is a problem, perhaps a freedof iterator can be created using similar logic to this PR? Potentially storing a Vector{UnitRange{Int}} for efficiency?
It is quite nice to have an easy way of doing e.g. r[fdofs] (even if not needed as it is used in the examples)

[KristofferC]

It is quite nice to have an easy way of doing e.g. r[fdofs] (even if not needed as it is used in the examples)

I guess I imagined something like https://github.com/JuliaData/InvertedIndices.jl where you would do r[Not(cdofs)].

[KnutAM]

Looks very neat - should for sure give the required functionality!
Nice to add to an example and docs if someone feels like removing ch.free_dofs!