New research on real-time deformable body dynamics (FEM).

[XperienS]

Download demos (zip, 1.1 Mb)
( controls: WASD + mouse [camera], F/G [apply force], Q [start simulation] )
This demos illustrates our research on deformable body dynamics.
We (me and FD) developed a method of simulating elastic and plastic deformation of bodies that can be easily integrated into a typical modern game rigid body dynamics pipeline.
The soft bodies are represented by particles and joints between them that can be solved together with rigid-body constraints.
This makes it possible to create constraints not only between the soft bodies, but also between soft and rigid bodies that will be solved simultaneously as any other constraint in the rigid body dynamics simulation.
As soon as your impulse-based dynamics engine is capable of handling points (that are totally the same as bodies, but without a rotational part in velocities, positions, forces and accelerations), and custom joints are supported, you can easily employ our method.
This approach finally makes it possible to generate contact points between soft and rigid bodies the same way they are generated for rigid bodies, making the interactions between soft and rigid bodies smooth, with the soft body being an integral part of the simulation.
The algorithm behind the soft body simulation is not a spring-mass simplification, but a true linear Finite Elements Method. The method uses tetrahedral finite elements, the infinitesmall strain theory with stiffness warping used to eliminate the artifacts of linearization, and the matrix form of the Hooke's law that enables the user to provide material properties in the form of the Poisson's ration and Young's modulus.
The algorithm is designed for real-time simulations; almost everything is precomputed and the code that is run every frame is as small as possible and easily parallelized per joint (per tetrahedron). Specifically there's no integration of a stiffness matrix or global matrix reassembling in this code.
For stiff systems, the accuracy of the typical Projected Gauss-Seidel solver is insufficient, so we propose a new method for solving LCPs that is based on conjugate gradients and has never been used in the game industry before; in fact the method was adapted by us so that with our modifications it is possible to substitute the most popular Projected Gauss-Seidel for this method in a typical impulse-based dynamics engine.
The new solver is also fit easily for GPU parallelization (as well as our FEM implementation), so we present you a version that runs on the GPU using CUDA technology.
Please try out the demos; if you are interested in the technology, feel free to contact us.

[fishboy82]

Seems good ,Could you write paper describe you new techinique?

[XperienS]

fishboy82
Thank you!
About the paper: right now we have only our old paper in Russian, and we will write a new one, in English (with additions and enhancements), as soon as possible. Although, this will require some time which we almost don't have now

[fishboy82]

That's great ,Look forward to see your paper,but also I am a little intersted to know the speed efficiency and memory need when solve the system

[XperienS]

Well, if I understood you correctly, here is data that you need:
1) solver complexity: O(n), n - number of constrained DoFs;
2) solver convergence: R-linear;
3) memory statistics:
Memory usage for body/node storage:
*per triple:
28 floats + 1 bool;
Memory usage for solver:
*per constrained DoF:
42 floats + 1 uint for PGS solver,
47 floats + 42 uints for our new solver;
*per triple:
6 floats for PGS solver,
6 floats + 2 uints for our new solver;
Memory usage per FE (tetrahedron):
111 floats;

The facts, needed to be considered:
* one rigid body is represented with two triples (linear and angular parts), while FEM node with only one (linear only);
* FEM tetrahedral constraint has 6 CDoFs and its own data to construct Jacobian (mem usage per FE);
* our new solver needs to know spectral radius of a matrix A, we used 10 iterations of Power Iterations method which is also O(n);

[fishboy82]

XperienS:
Thanks alot for you reply,that's the data I need, Actually I am more intersting about your new CG LCP sover, so you mean you solve rigid body and sift body all use your new solver? Or you use traditional PGS solver for rigid body and use your new solver for rigid-soft constraint(contact etc)？

[XperienS]

fishboy82
We're using one solver for solving the whole system, either PGS or our new CG-based one. We're not combining them, that's the point of our research - to make one system which could be solved at once, w/o any tricky coupling techniques. So yes, all the constraints (between rigid body, rigid-soft body, and between soft bodies) are solved within one system with a chosen solver.
Data for PGS in my previous message was presented for comparison of memory usage to our new solver.

[Dirk Gregorius]

One nice thing about PGS is that with each iteration the result becomes better. With CG can you get a worse result when you go the next iteration. You basically get spikes in the convergence graph. Another nice thing about PGS is that it can be warmstarted so you get an amortized iteration count. Usually games use 5-10 iterations. You seem to need iterations in the range of ~50.

Can your solver be warmstarted?
Does your solver converge gradually?
What is the typical iteration count and how would you compare it to one pure PGS iteration?

I have seen combinations of PGS and CG. For example Moten Silcowitz from Diku (author Jinnengine http://code.google.com/p/jinngine/ ) uses a Fletcher–Reeves type nonlinear non-smooth conjugate gradient (NNCG) type method. How does your method compares to his research?

Cheers,
-Dirk

[fishboy82]

Hi Dirk:
have you got the pdf of the paper<<A nonsmooth nonlinear conjugate gradient method for interactive contact force problems >>,If you have, could you send me one?

[Dirk Gregorius]

Please ask Morten for the paper. He is very responsive.

[XperienS]

Dirk Gregorius
> Can your solver be warmstarted?
Yes, it definitely can be simply warmstarted (as well as CG though) and it amortizes iteration count too.

> Does your solver converge gradually?
Our solver have the similar behaviour here as CG solver. It is not 'unpredictable' - it makes steady progress towards the solution, although there could be 'spikes' in the Euclidean norm of the residual that you mentioned.

> What is the typical iteration count and how would you compare it to one pure PGS iteration?
Actually, 'typical iteration count' varies widely: e. g. 3 swinging chains of 50 ball-in-socket joints (right end of its initial position fixed) could converge simply within 15 iterations; but our FEM tests usually require 30-150 iterations (150 when stressed hardly), while PGS could easily introduce instabilities even within 300-1000 iterations (when used on highly stiff systems). But one iteration of CG-based solver is 3-4 times heavier than one PGS iteration, although CG-based solver is highly parallel.

> How does your method compares to his research?
We haven't read his research yet, we'll try to make a comparison later in a spare time; but in contrast to what you've said, our approach has nothing in common with PGS.

[Dirk Gregorius]

Thanks for the information. I would be interested to read about your research .

...while PGS could easily introduce instabilities even within 300-1000 iterations (when used on highly stiff systems)...

This is not due to PGS, but due to the stabilization method when solving on the velocity level.

[nisse]

Impressive work! Are you solving the contacts as unilateral contacts - with friction - using the CG solver as well? Do you have your own method for handling the active set in the LCP solve or do you use some standard method? That is the hardest part to parallelize as I understand the problem.

Keep up the good work - and please publish when you feel ready for it

[fishboy82]

Hi Dirk：
I read the paper <<A nonsmooth nonlinear conjugate gradient method for interactive contact force problems >>
But there seems some problem in the paper . As they say r = (H − I)Lamda+h can be
interpreted as the gradient of a non-smooth nonlinear quasi-quadratic function f(Lambda) ≡ 1/2|| r||^2 and the gradient of this function is gradf = -r;but when I deduce this gradient actually should be (H-I)*r rather than -r have you deduced that?

Nisse and XperienS：
I only know to use CG to solve linear eqaution but don't know how to use it solve LCP,
Could you refer me some Material how to use active set in CG method to solve LCP with unilateral contacts and box constraint like friction？

[spade]

Hi XperienS,
is your method able to solve softbody (FEM) and contact constraints (not only joints) between the soft and rigid body together in a coupled system?
In my engine I've tried to solve it and I've end up with a linear system (LCP) in which the matrix is still symmetric but no more positive defined.
Indefinite systems shouldn't be solved by a simple CG method but require other algorithms and I believe that PCG doesn't work for a indefinite LCP.
Have you encountered this problem?


I'm very interested on your research paper