Ageia paper on position based physics

KenB
Posts: 49
Joined: Sun Dec 03, 2006 12:40 am

Regularized methods

Post by KenB »

Someone might be interested in looking at this:
http://www.vrlab.umu.se/research/flexible.shtml

And specifically at:
Interactive simulation of elastic deformable materials at
http://www.vrlab.umu.se/research/flexible/ISEDM.pdf

This implementation has not been optimized for speed, but mainly for evaluating the method, constitutive equations etc. The solver can be optimized way beyond the current implementation.

/Kenneth
bone
Posts: 231
Joined: Tue Feb 20, 2007 4:56 pm

Re: Regularized methods

Post by bone »

KenB wrote: And specifically at:
Interactive simulation of elastic deformable materials at
http://www.vrlab.umu.se/research/flexible/ISEDM.pdf
I looked at it. As they note, it's really two parts: one is a more physically accurate way to model elastic tetrahedrons, the other is a new integration scheme.

I'm still trying to fully understand exactly what they're doing in the new integration scheme and wondering if it can be applied successfully in a more general case, for example an iterative solver.
jiangwei
Posts: 23
Joined: Wed Nov 30, 2005 11:07 am
Location: China

Post by jiangwei »

Is there any open source physics engine which used position based method?
thanks
raigan2
Posts: 197
Joined: Sat Aug 19, 2006 11:52 pm

Post by raigan2 »

jiangwei wrote:Is there any open source physics engine which used position based method?
thanks
I'm working on one in 2D, however I've recently run into a problem: http://www.continuousphysics.com/Bullet ... php?t=1152

The solver is excellent -- super-stable just like Jakobsen's, but much more applicable as you can extend it to rigid bodies, motors, etc..
JohnJensen
Posts: 2
Joined: Tue Oct 02, 2007 12:54 am

Re: Ageia paper on position based physics

Post by JohnJensen »

Hello-

I'm interested in implementing this paper for a cloth simulation I'm working on. The problem is that the "solver" section is fairly vague. Well, it's not vague, but it implies knowledge about Gauss-Seidel which -- even after a bit of research -- I cannot find a clear example of its usage. I studied numerical methods a bit in college, but I'm not familiar with all the intimate details of course.

I'm confused as to whether or not I'm supposed to be creating a giant matrix and solving that, or if Gauss-Seidel is just a nice name for a relatively simple practice (with a strong mathematical concept "behind the scenes").

I'm mostly looking to implement the cloth solution from this paper. Any direction, code samples, etc. should be useful. If someone were to walk me through an iteration of a basic 2x2 piece of cloth, I'd be very thankful

Code: Select all

p2 -- p4
| \    |
|    \ |
p3 -- p1
So say I have edge constraints:

C1(p2, p3), C2(p3, p1), C3(p1, p4) C4(p4, p2)

Should I also have diagonal constraints? (e.g. C5(p1, p2) & C6(p4, p3)) Reading the ppaer, I believe I should have C5(p1, p2) but not the other diagonal since it's not a triangle edge. Am I reading that correctly? Why not have the other diagonal?

And I know I'll have one bend constraints Cbend(p1, p2, p3, p4) between two two triangles /\(p1, p3, p2), /\(p1, p2, p4).

So, I have those 5-7 constraints.. How do I solve this thing?

Am I reading too heavily into this? I have a feeling I am.

Am I supposed to just go from constraint C1 through CN and "apply" (project) the constraint to p1-p4.

So C1 updates p2 & p3, C2 updates p3 & p1 (using the results from C1), and so on? In the case of cloth, this would be extremely similar to the Jacobson method except adding bend and collision constraints?

Is the real novelty here for cloth is that we are using predicted positions to modify the velocity ala the Nonconvex Rigid Bodies w/ Stacking paper (albeit in a non-impulse way)?

I'm new to a lot of this, so any help or guidance would be appreciated. While I value the mathematical theory behind things, all I've been able to find is theory so some concreteness would help!

Thanks,
-john
raigan2
Posts: 197
Joined: Sat Aug 19, 2006 11:52 pm

Re: Ageia paper on position based physics

Post by raigan2 »

JohnJensen wrote: Am I supposed to just go from constraint C1 through CN and "apply" (project) the constraint to p1-p4.

So C1 updates p2 & p3, C2 updates p3 & p1 (using the results from C1), and so on? In the case of cloth, this would be extremely similar to the Jacobson method except adding bend and collision constraints?

Is the real novelty here for cloth is that we are using predicted positions to modify the velocity ala the Nonconvex Rigid Bodies w/ Stacking paper (albeit in a non-impulse way)?
You're correct -- you just iterate over all constraints, projecting one-by-one. The way that positions are updated (using "current" and "predicted") are in fact identical to Jakobsen ("oldpos" and "newpos") if you work it through on paper.

The main benefit of this paper is that it shows you how to properly extend Jakobsen's method to any sort of constraint (provided it's described as a function of particle positions). Personally I struggled for a long time with getting angle constraints to be stable, and the solver presented by Muller works excellently for this (in 2d at any rate).
JohnJensen
Posts: 2
Joined: Tue Oct 02, 2007 12:54 am

Re: Ageia paper on position based physics

Post by JohnJensen »

Thanks for your help.. It's getting clearer by the minute.

So is the "solver" the section where he described constraint projection that conserves momentum? E.g, is the conservation the "condition" he needs to derive a numerically stable solver?

So what's the "Gauss-Seidel-type iteration"? Is it "merely" the process of repeatedly applying the constraints?

What about diagonal constraints, are there one or two in the sample I presented above?

Sorry I'm so questiony-questiony -- I am experimenting with this at the moment, I just want to make sure I'm understanding things.

Thanks!
Antonio Martini
Posts: 126
Joined: Wed Jul 27, 2005 10:28 am
Location: SCEE London

Re: Ageia paper on position based physics

Post by Antonio Martini »

raigan2 wrote:The main benefit of this paper is that it shows you how to properly extend Jakobsen's method to any sort of constraint (provided it's described as a function of particle positions). Personally I struggled for a long time with getting angle constraints to be stable, and the solver presented by Muller works excellently for this (in 2d at any rate).
a couple of alternative methods for solving the problem in a more general way are mentioned here:

http://www.bulletphysics.com/Bullet/php ... f=4&t=1593

cheers,
Antonio
raigan2
Posts: 197
Joined: Sat Aug 19, 2006 11:52 pm

Re: Ageia paper on position based physics

Post by raigan2 »

Thanks -- i've been reading the thread (and papers) with great interest. Sadly I've only taken one year of linear algebra so I'm not entirely confident about my understanding..

I'm especially interested in the equivalent row-by-row/SCAAT version (as opposed to the big-matrix approach) since it seems like it would end up looking almost the same as the PBD solver except that it would handle masses/weighting differently; currently I'm using an invented scaling term since the revised paper's formula led to very slow convergence.

I've been meaning to post about my simulator, as increasingly the main problem has been collision detection (or, updating jacobians/attach points per iteration). My current solution is very similar to the Stewart/Trinkle or Anitescu/Potra solvers that have been posted here from time to time (as least, as far as I understand) -- unfortunately this results in a large number of constraints being generated.
Antonio Martini
Posts: 126
Joined: Wed Jul 27, 2005 10:28 am
Location: SCEE London

Re: Ageia paper on position based physics

Post by Antonio Martini »

raigan2 wrote:Thanks -- i've been reading the thread (and papers) with great interest. Sadly I've only taken one year of linear algebra so I'm not entirely confident about my understanding..

I'm especially interested in the equivalent row-by-row/SCAAT version (as opposed to the big-matrix approach) since it seems like it would end up looking almost the same as the PBD solver except that it would handle masses/weighting differently; currently I'm using an invented scaling term since the revised paper's formula led to very slow convergence.
you can try to parametrize the convergence rate by redefining the mass weights the following way.

w1 = (m2+q)/(m1+m2+q+q)
w2 = (m1+q)/(m1+m2+q+q)

for q = 0 you get the usual weighting
for q = large number you get approximately 0.5 for both weights
q can be a global parameter.

the idea is to start with a low q and increase it at each iteration, so you get the mass difference effect and a faster convergence. This is what the Q covariance matrix in the Kalman filter does. So referring to my other post(i have updated it slightly) this is the same as not having equation (e) and just updating the inverse mass at every frame by eq(b). By including equation (e) the convergence is even faster. The total error may still be the same, however the maximum local error should be smaller depending on the choice of Q. So we dont get localised overstretched edges. My tests until now have been purely at numerical level(comparing algorithms by running them on the same data) with no graphical output however the weighting above make sense also at intuitive level.

cheers,
Antonio
Last edited by Antonio Martini on Thu Oct 25, 2007 12:11 am, edited 1 time in total.
raigan2
Posts: 197
Joined: Sat Aug 19, 2006 11:52 pm

Re: Ageia paper on position based physics

Post by raigan2 »

That's a really interesting idea!

I was referring to the application of the weights -- using them to scale the projections. My (possibly mistaken) understanding of this process is that the weights determine how far we descend along each particle's gradient, acting like the damping parameter in a damped-Gauss-Newton optimization/least-squares-fitting process.

The individual particle weights are as in the paper: w_i = 1/m_i

We project each particle along its gradient by a certain amount: dp_i = -step_i*dC_i
where dp_i is the projection to be applied to particle i, and dC_i is the gradient of the constraint function wrt the particle.

The original paper suggests this formula for the step length (apologies for uglyness):
step_i = ((n*w_i)/Sum_j(w_j)) * (C / Sum_j( Dot( dC_j, dC_j) ))

Where C is the current value of the constraint function, and n is the "cardinality" of the constraint (the number of particles). This formula doesn't work, it leads to overshooting/oscillation when the mass ratios aren't 1:1.

The revised PBD paper gives the formula as:
step_i = w_i * (C / Sum_j( w_j * Dot( dC_j, dC_j) ))

I've tried implementing this three times, and the results have always been the same: very slow convergence.
I tried using the following variation, however it didn't work (things exploded):
step_i = w_i * (C / Sum_j( Dot( w_j*dC_j, w_j*dC_j) ))

The formula that I'm currently using is the result of a mistake, when implementing the original paper I tried scaling the gradients _before_ calculating their squared length: let s_i = (n*w_i)/Sum_j(w_j)
step_i = s_i * (C / Sum_j( Dot( s_j*dC_j, s_j*dC_j) ))

This formula provides the best scaling so far; it never overshoots (even with mass ratios of 1:100), and provides very fast convergence.

I'm sure that there is a more optimal/correct scaling factor, since the use of the cardinality term is totally ad-hoc and shouldn't be necessary. For instance, you should be able to consider each component of the particle position seperately (i.e convert the constraint function C(p1,p1) to one which takes scalars C(p1x,p1y,p2x,p2y)) and the results should be identical, however since there are 4 scalars rather than 2 vectors, this changes the cardinality and the resulting projections.. that doesn't seem right to me.

I may be incorrect as far as construing the PBD solver to be doing a least-squares-fit/optimization process, the formulas simply seem very similar to Newton-type optimization. I would _really_ like to understand why my formula works so much better, and what the optimal formula is -- this is what interests me most about your Kalman filter method, it seems like a more correct method (no magical "cardinality" values floating around). I haven't wrapped my head around it yet, unfortunately "filter" seems very similar to "kernel" in that it seems to be only an abstract engineering tool and not something that can be easily visualized!
Antonio Martini
Posts: 126
Joined: Wed Jul 27, 2005 10:28 am
Location: SCEE London

Re: Ageia paper on position based physics

Post by Antonio Martini »

raigan2 wrote: I may be incorrect as far as construing the PBD solver to be doing a least-squares-fit/optimization process, the formulas simply seem very similar to Newton-type optimization. I would _really_ like to understand why my formula works so much better, and what the optimal formula is -- this is what interests me most about your Kalman filter method, it seems like a more correct method (no magical "cardinality" values floating around). I haven't wrapped my head around it yet, unfortunately "filter" seems very similar to "kernel" in that it seems to be only an abstract engineering tool and not something that can be easily visualized!
equations (a-e) of my post dont require any filter notion, you can just plug in the values and and even keep them in matrix form. If you get an improvement then they can be written in a much simpler form as many terms typically cancel out.

for my tests i usually use newmat:

http://www.robertnz.net/nm_intro.htm

it's a very simple to use and flexible library. However it's terribly slow, when im happy with the results then i simplify the computations and move on other far more efficient libraries.

if you want to try the EKF approach just send me a private message so we can get it to work.

cheers,
Antonio
Last edited by Antonio Martini on Wed Oct 24, 2007 11:48 pm, edited 3 times in total.
raigan2
Posts: 197
Joined: Sat Aug 19, 2006 11:52 pm

Re: Ageia paper on position based physics

Post by raigan2 »

AntonioMartini wrote: equations (a-e) of my post dont require any filter notion, you can just plug in the values and and even keep them in matrix form. If you get an improvement then they can be written in a much simpler form as many terms typically cancel out.
Thanks very much, I'm planning on working through this tomorrow, hopefully I'll be able to manage. I feel quite stupid compared to all of the experts on this board, I'm definitely more of a "technical artist" than a proper programmer in that my math background is weak.

Dirk Gregorius wrote:Strange that you have these issues. The formulas seem correct. I derived the formulas myself and in the new version it correct now. The derivation is extremely simple (using block vector notation):
I asked Mr. Muller about this issue, and he said the revised version is working for him, and that my formula might ismply be introducing a successive-overrelaxation-type of bias. However, I've received a private message from someone who was also having slow convergence with the revised formulas.

The revised formula may very well be correct; the issue may be with my expectation of how it should perform -- my formula might simply be more aggressively descending the gradient, and it turns out that for the types of constraints I'm using this bias never results in problems (overshooting/etc), thus I get accelerated convergence. Currently, I get very solid behaviour with a low number of iterations (4 iterations for this biped: http://www.metanetsoftware.com/robotolo ... ktest.html); you can violently swing it around and the bodies stay together -- even two iterations is enough to prevent obvious visual gaps.

However with the revised formulas, 4 iterations feels very weak/loose -- if you move one of the bodies across the screen, after the solve the remaining bodies will be spread out across the screen.. huge gaps. At around 10-20 iterations it becomes solid -- perhaps this is normal behaviour, I know box2D uses 10s of iterations. So perhaps the only problem is my high expectations of convergence; nevertheless, the fact that my current method behaves much better is something I'd like to understand! Also, since with the PBD approach you have to rebuild the gradients/jacobians at each iteration, keeping the iteration count low is desirable.

There is also the possibility that my code is flawed; once I have time to return to this I'll clean up and post some of the code.
raigan2
Posts: 197
Joined: Sat Aug 19, 2006 11:52 pm

Re: Ageia paper on position based physics

Post by raigan2 »

Working through a-d (skipping e, and using Q=0, R=0) I get eq9 in the PBD paper; this is the same formula as Dirk's most recent post.

If I understand things correctly, my ad-hoc formula is equivalent to using a non-zero value for R; I'm going to try a sanity check implementation of a simple chain of jointed rigid bodies.
Antonio Martini
Posts: 126
Joined: Wed Jul 27, 2005 10:28 am
Location: SCEE London

Re: Ageia paper on position based physics

Post by Antonio Martini »

raigan2 wrote:Working through a-d (skipping e, and using Q=0, R=0) I get eq9 in the PBD paper; this is the same formula as Dirk's most recent post.

If I understand things correctly, my ad-hoc formula is equivalent to using a non-zero value for R; I'm going to try a sanity check implementation of a simple chain of jointed rigid bodies.
if they are equivalent you should be able get almost identical numerical results. If thats the case you may want to test the effect of Q and (e) on the final result.

cheers,
Antonio
Post Reply