Alternative to verlet integration

Started by
7 comments, last by Alessandro 12 years ago
I'm developing an application that allows to model hair: each strand is made of 17 vertices and I use verlet integration to constraint vertices so that when vertices are moved the distance between each vertex remains the same.
It works just fine, but it gets slow when there are many hair selected. That's caused by the sqrt() needed in the verlet function.

So I'd like to ask if there might be some alternatives in order to have a faster solution. Thanks
Advertisement
I think you are talking about constraints, not verlet integration per se (although they are often used together).

I don't think you can save yourself the sqrt, but you can probably get away with some dirty approximation. See if you can use a polynomial approximation or if you can use this type of trick.
It's probably not the verlet integration (Jacobsen style, I presume?) that's causing the slowdown. It's all the square root calls. You could avoid those all together by defining your hair strains as series of 1d springs instead (one for each dimension). The force applied by a 1d damped spring is:

force = - ( stiffnes * (length - rest length) + damping * velocity )

Since this is 1d you don't need pythagoras theorem - and hence no square root - to find the spring's length. No need for a normalized vector or vector projection of velocity for the damping part either. They are all scalar values. By defining unique rest lengths for each particle pair for each dimension, you get a hair strain with built-in linear and angular springs. Feel free to ask for more details if my explanation is a bit unclear.

Also, by using a Catmull-Rom cubic spline function to render your hair simulation, you can reduce the number of vertices by more than half and still keep the soft, rounded look of the strains.

Cheers,
Mike

I think you are talking about constraints, not verlet integration per se (although they are often used together).

I don't think you can save yourself the sqrt, but you can probably get away with some dirty approximation. See if you can use a polynomial approximation or if you can use this type of trick.


Hi Alvaro, thanks for the prompt answer. Well, I apologize because I didn't specify that I was already using that "babylonian method" (I think is called like that).
<p>&nbsp;</p>
<p>
<br />
It's probably not the verlet integration (Jacobsen style, I presume?) that's causing the slowdown. It's all the square root calls. You could avoid those all together by defining your hair strains as series of 1d springs instead (one for each dimension). The force applied by a 1d damped spring is:<br />
<br />
force = - ( stiffnes * (length - rest length) + damping * velocity )<br />
<br />
Since this is 1d you don't need pythagoras theorem - and hence no square root - to find the spring's length. No need for a normalized vector or vector projection of velocity for the damping part either. They are all scalar values. By defining unique rest lengths for each particle pair for each dimension, you get a hair strain with built-in linear and angular springs. Feel free to ask for more details if my explanation is a bit unclear.<br />
<br />
Also, by using a Catmull-Rom cubic spline function to render your hair simulation, you can reduce the number of vertices by more than half and still keep the soft, rounded look of the strains.<br />
<br />
Cheers,<br />
Mike<br />
</p>

Hi Mike, thanks for your detailed answer. Well, I do understand the concept (working on 1D strand doesn't imply using square roots) but I'm not sure how'd I modify my code. It would great if you could give me a hint about how to modify it in order to use 1D springs.
Here is a excerpt from my verlet routine (and yes, verlet integration was a wrong term):


constraint_distance=myGuideHair.hair_constraint_distance; // the constraint distance between each point

//This Section Moves All Spline Points Except Points ( 0 , 1 And Last Point )
//Get Two Points In Spline And Create A Vector Between Them
//Calculate The Difference Between The Actual Distance And The Distance They Should Be
//Divide By Two Because We Are Moving Two Points Either Towards Or Away From Each Other.
//So We Move Each Point Only Half The Distance. But In Opposite Directions
for (int p=1; p<HAIRPOINTS-1; p++)
{
verletVec.x = myGuideHair.hairControlPoint

.x - myGuideHair.hairControlPoint[p+1].x;
verletVec.y = myGuideHair.hairControlPoint

.y - myGuideHair.hairControlPoint[p+1].y;
verletVec.z = myGuideHair.hairControlPoint

.z - myGuideHair.hairControlPoint[p+1].z;

magnitude = fastSqrt_Bab_2(verletVec.x*verletVec.x + verletVec.y*verletVec.y + verletVec.z*verletVec.z); //Calculate The Distance Between Both Points

difference = (magnitude - constraint_distance) / 2.0;
//Normalize The Vector To A Unit Vector
if (magnitude==0) magnitude = 1;
verletVec.x = verletVec.x / magnitude;
verletVec.y = verletVec.y / magnitude;
verletVec.z = verletVec.z / magnitude;
//Increase The Length Of The Vector By The Difference ( This Distance Is Used To Move Both Points To wards Or Away From Each Other )
verletVec.x*=difference;
verletVec.y*=difference;
verletVec.z*=difference;
//Move This Point In A Direction
myGuideHair.hairControlPoint

.x-=verletVec.x;
myGuideHair.hairControlPoint

.y-=verletVec.y;
myGuideHair.hairControlPoint

.z-=verletVec.z;
//Move This Point In The Opposite Direction To The Previous
myGuideHair.hairControlPoint[p+1].x+=verletVec.x;
myGuideHair.hairControlPoint[p+1].y+=verletVec.y;
myGuideHair.hairControlPoint[p+1].z+=verletVec.z;
}



I should also mention that each hair, made of 17 points, is done so that point 1 is actually the "hook point" (i.e. the point attached to the head mesh; point 0 and last point are the first and last control point of the catmull rom spline...
FYI for other posters, I think this thread might give some helpful background information on how he's doing this (I'm pretty sure he's using the same method I used in the sample code).

I don't know of a way to do this without the square root, so I'm sorry I can't suggest a root-less alternative. However, if you wanted to dig into SIMD operations, like SSE, and if you write the instructions in an efficient manner, you can make it run tons faster. If this is a CPU intensive program, I'd seriously consider optimizing bottlenecks with SIMD code.
[size=2][ I was ninja'd 71 times before I stopped counting a long time ago ] [ f.k.a. MikeTacular ] [ My Blog ] [ SWFer: Gaplessly looped MP3s in your Flash games ]

FYI for other posters, I think this thread might give some helpful background information on how he's doing this (I'm pretty sure he's using the same method I used in the sample code).

I don't know of a way to do this without the square root, so I'm sorry I can't suggest a root-less alternative. However, if you wanted to dig into SIMD operations, like SSE, and if you write the instructions in an efficient manner, you can make it run tons faster. If this is a CPU intensive program, I'd seriously consider optimizing bottlenecks with SIMD code.


Hi Cornstalk that's correct, the code I'm using uses the same method of the one that you kindly provided and that was very helpful. I really know nothing about SIMD code, I assume is something like assembly.
Hi Alessandro,

I'm afraid I can't simply change a few things in your code snippet in order to make it work like I suggested. This is what I recommend you to do: You already have a hair strain class called myGuideHair and a hair strain "particle" class called hairControlPoint. I suggest you create a new class that manages the constraint between two control points. This class contains a pointer to each of the two control points and a unique rest length vector containing the three 1d rest lengths defining the length and angle of this particular section of hair strain. The class also contains a method for applying forces or impulses to the two control points in order to make them stay the predefined distances apart. Exactly what method you use - Thomas Jacobsen "verlet" constraint, Hooke's law damped spring etc. - is less important. for each dimension of the hair strain section calculate the displacement = control point distance - control point rest distance. Apply a proportional force or impulse in that dimension only.
The only drawback with this method is that you cannot rotate the hair strains. If you need to rotate them, you will need to define the 1d rest lengths in some local space and then rotate that space. I have made an example of how to implement this in c++ that you may download and use. In my example, I am using the springs for a soft body implementation with pressure, so simply ignore this part and look at the springs only.

Cheers,
Mike
Hi Mike, thank you very much for the explanation, and all the time you spent for it and to make the example that you kindly shared. I'm going to study it and see if I can make it!

This topic is closed to new replies.

Advertisement