[SofaCUDA] ElementFEMForceField: Generic CUDA implementation#6071
[SofaCUDA] ElementFEMForceField: Generic CUDA implementation#6071fredroy wants to merge 21 commits into
Conversation
fredroy
added
pr: enhancement
alxbilger
left a comment
alxbilger
left a comment
There was a problem hiding this comment.
A criticism is that the code supposes that the space is in 3D. There are also some dynamic checks on the number of nodes per element. But this data is constexpr.
fredroy
force-pushed
the
femelasticity_cuda
branch
from
4f653db to
e41d0b6
Compare
Taken into account by 🤖 |
Extract assembled stiffness matrices into a separate contiguous buffer (m_assembledStiffnessMatrices) to replace getReadAccessor calls on Data<vector<FactorizedElementStiffness>> inside parallel forEachRange lambdas. The read accessor acquires a shared lock on the Data object, causing contention across threads and effectively serializing the parallel work during CG iterations. Using a direct const reference to a plain vector eliminates this synchronization bottleneck (~3x speedup in parallel mode). As a secondary benefit, the contiguous buffer only stores the assembled 24x24 matrices (~4.6 KB each) rather than the full FactorizedElementStiffness structs (~14 KB each), improving cache utilization.
Replace hardcoded 3D assumption and extern "C" + switch(nbNodesPerElem) runtime dispatch with fully compile-time C++ template parameters <T, NNodes, Dim>. All kernel dimensions, stiffness block sizes, and gather loops are now generic over Dim. The .inl callers use a single template call with constexpr nNodes and dim from the trait, eliminating both the if-constexpr type branching and the runtime NNodes switch. Explicit template instantiations in the .cu files provide the needed symbols. Applied to both ElementLinearSmallStrainFEMForceField and ElementCorotationalFEMForceField CUDA implementations.
This reverts commit e41d0b6.
fredroy
force-pushed
the
femelasticity_cuda
branch
from
6a2d395 to
b6c5e82
Compare
|
@courtecuisse this is a PR on which you could make a review since these FEM models are meant to replace old FEM constitutive laws later on. |
|
So I've tested it. I've actually tested it in float, I don't know why I went down this road but hey, if it passes in float then it is also good in double. As a recall, I wanted to validate this implementation by using regression tests. The idea was to generate the CPU regression reference and then use it as reference for the CUDA execution and see what is the result. So I selected 4 scenes to test is against the CPU version :
To make those scenes run properly in float, I had to scale them by 1e3 to avoid precision breakdown.
Which is better than the threshold for the CPU tests. And recall that those scenes are scaled by a factor 1e3, meaning that, those precision values should actually be around 1e-10 (if double where used). So for me this is good. |
4 participants
Based on
Implementation of the #5882 for SofaCUDA.
Claude did the heavy-lifting job :
And the human implemented examples, comparison and benchmarks with the current implementation of FEM in SofaCUDA.
As for the bench (only for tet and hex):
In a nutshell,
Benches (only corotational...)
(i7 13700K + 4080Ti)
can be launched like that :
LinearSolver is CG, 250 steps, tol=1e-12
40x10x10 grid (1000 steps, 4000 nodes, 3159 hexa )
76x16x16 grid (1000 steps, 19456 nodes, 16875 hexa ) (no cpu)
For the fun, the hexa new version with the a 200x30x30 (180k points, 167k hexa) and a downgraded CG (50 steps, tol of 1e-06) :
1000 iterations done in 65.5718 s ( 15.2504 FPS)By submitting this pull request, I acknowledge that
I have read, understand, and agree SOFA Developer Certificate of Origin (DCO).
Reviewers will merge this pull-request only if