For 2d brittle-fracture problems, TNNMG is too slow. For load steps away from the actual fracture, the operator-splitting algorithm is twice as fast, if not faster.
I measured the TNNMG performance with the gprof profiler. About a quarter of the time is spent in each of the three methods: the presmoother,
IntegralFunctionalConstrainedLinearization::bind. I think there a still ways to make the code quite a bit faster; here are patches and ideas:
IntegralFunctionalConstrainedLinearization::addRowTransposedXMatrixXRow: Exploit that the result matrix is symmetric.
IntegralDirectionalRestriction::subDifferential: Precompute the product
evaluationMatrix * direction.
Use gradient-aware truncation. Early tests suggest that this reduces the number of iterations at least when the load is low.
Use a different line search algorithm. Calls to
subDifferentialare not cheap (though much less expensive with the patch mentioned above), and we know that the functional is smooth.
Exploit matrix symmetry in
evaluationMatrixto contain derivatives of scalar shape functions, rather than lineared strains. I believe that several parts of the algorithm are limited by the memory speed, rather than by computing speed. Storing only shape function gradients instead of strains would reduce the required memory bandwidth considerably. It would also reduce our overall memory requirements which, according to Daniel's measurements, are nothing to brag with either.
Even the entries of the
evaluationMatrixare sparse. In 2d, the entries are 4x3 matrices, but 7 of the 12 are always zero. If the matrix-vector multiplication becomes a bottle-neck then using a dedicated matrix type may speed up the code.
SymmetricMatrixobjects for small symmetric objects. Not sure how much this would buy us.