sofa-framework / sofa

Real-time multi-physics simulation with an emphasis on medical simulation.
https://www.sofa-framework.org
GNU Lesser General Public License v2.1
922 stars 311 forks source link

MechanicalMatrixMapper fast matrix product does not detect changes in the jacobian matrix #2687

Open alxbilger opened 2 years ago

alxbilger commented 2 years ago

In MechanicalMatrixMapper, the fast matrix product relies on the fact that both sparse matrices do not change their structure. If one of the matrices structure changes, the fast matrix product is no longer valid and leads to wrong results. One case where the matrix structure changes is with topological changes. This is supported in MechanicalMatrixMapper because there is a check that the topology changed compared to the previous time step. In that case, the matrix intersections is invalidated and it triggers a new computation of the intersection. Other cases are not supported. For example, @Camille-K experienced the issue: the mapping jacobian matrix structure changed because of a user interaction, moving a frame from a beam to another. More generally, it could be profitable to have a generic mechanism detecting that a matrix structure changed.

ljluestc commented 9 months ago

class MatrixStructure {
public:
    MatrixStructure() : version(0) {}

    void updateStructure() {
        // Update the matrix structure (e.g., due to user interaction)
        // You can update indices, dimensions, or any relevant information here
        version++;
    }

    int getVersion() const {
        return version;
    }

private:
    int version;
};

// Usage example:
MatrixStructure matrixStructure;

// Update the matrix structure when changes occur
matrixStructure.updateStructure();

// Check if the matrix structure has changed before performing a computation
int currentVersion = matrixStructure.getVersion();
if (currentVersion != previousVersion) {
    // The matrix structure has changed, perform necessary updates
    // For example, recompute matrix intersections or mappings
    // ...

    // Update the previous version to the current one
    previousVersion = currentVersion;
}