# Deformation matrix, Math hack for stability on large simulation steps?

So there is a numeric technique for updating a deformation gradient in MPM that goes as:

$$F_{n+1} = (I + \nabla \vec v \Delta t)F_n$$

This works for small time steps but for large time steps numerical issues ensue. One of the biggest issues is that I assume that the deformation matrix always has a determinant of one, but for large time steps it can flip the sign.

I am looking for a numerically stable (accuracy does not matter) way to modify the matrix. Basically if the determinant is negative, I want to tweak the resulting $$F_{n+1}$$ matrix such that its determinant is now positive, and it remains stable in terms of the overall simulation.

The discrete formula comes from the differential equation:

$$\frac{dF}{dt} = \nabla \vec v F$$

I am not sure if I can leverage the analysis to self correct for the determinant.

• Try implicit Euler, you'll have to solve a linear system however. For very large time steps this may be cheaper than doing many small time steps with explicit Euler however. I have to idea what objects F and v are in your question by the way, from what I read $F$ seems to be a matrix, but what is $\nabla \vec{v}$ then? Is it some kind of Jacobian matrix? Commented Jul 7, 2023 at 18:42
• $F$ is the deformation gradient, $\nabla \vec v$ is the gradient of the velocity field (similar to a jacobian matrix). I cannot do implicit Euler, it takes far too much time. Commented Jul 7, 2023 at 20:03
• If you have multiple steps you can try to apply something like fast explicit diffusion: mia.uni-saarland.de/Research/SC_FED.shtml or the scheduled jacobi relaxation method. I don't know whether it will work well for your problem though. Beyond that you can make the resolution of your grid coarser (if your matrices arise from a grid discretisation). Is $\nabla v$ fixed by thevway or does it change between steps? What is it a function of? What did you use for the implicit system solution that it was so slow? Have you tried a Jacobi preconditioner instead of just Richardson? Commented Jul 7, 2023 at 21:50
• I am doing MLS-MPM, simulating a fluid with 50k particles doing one simulaiton step per frame at 60 fps. $\nabla v$ changes between time steps it;s the velocity of the particles when transported ot the grid. I did not implement it, I read an article by the creator of the MLS-MPM method which said it's likely not even worth trying. Commented Jul 7, 2023 at 22:04
• If $\nabla v$ is a matrix of the velocities of the particles how is it square? Isn't a velocity a $d$-dimensional vector and you have $N$ particles, so $d\times N$? Anyways assume for a second that $\nabla v$ was a fixed matrix $-A$, then the iteration that you have is convergent for $\|I-\Delta t A\|<1$, that's generally where the step size constraint comes from, for a symmetric positive definite $A$ and the $2$-norm this generally translates to $\Delta t \in (0, 2/\lambda_{max})$. You can't avoid this with explicit Euler, unless you consider a multistep relaxation method as mentioned. Commented Jul 7, 2023 at 22:34

## 1 Answer

The part about your 'stability' requirements seems impossible to answer without more detail. But if you just want to ensure that $$\det F_{n+1} = 1$$ at each step, then a quick fix is scaling your update rule to ensure that this property holds: $$F_{n+1} = \alpha(I + \nabla \vec v \Delta t)F_n,$$ where $$\alpha>0$$ is chosen so that $$\det \alpha(I + \nabla \vec v \Delta t) = 1.$$ If $$\nabla \vec v$$ and $$\Delta t$$ are vectors, you can compute the determinant quickly with the matrix determinant lemma. Then don't forget to scale correctly: for an $$m\times m$$ matrix $$A$$, $$\det \alpha A = \alpha^m \det A$$.

• if m is squared this means that if the determinant is negative then then there's no way to find $\alpha$ such that $\alpha A$ = 1. since ti will always be negative. Commented Jul 7, 2023 at 18:30
• True (I assume you mean if $m$ is even). I have no idea how your simulation works, but I assumed that the determinant got negative because of the combined uncorrected effect of several time steps, not just of a single one. Commented Jul 7, 2023 at 18:47