Notes 1D Elastic Continuum � From last class: elastic rod � Required reading: • linear density “rho” (not necessarily constant) • Baraff & Witkin, “Large steps in cloth • Young � s modulus E (not necessarily constant) animation”, SIGGRAPH � 98 • Paratemerized by p • Grinspun et al., “Discrete shells”, SCA � 03 • Bridson et al., “Simulation of clothing with � � � � � � x ( p ) = 1 � p x ( p ) � 1 �� folds and wrinkles”, SCA � 03 � p E ( p ) � � � � � � � � � � If homogenous, simplifies to: � 2 x � 2 x � t 2 = E � � p 2 cs533d-term1-2005 1 cs533d-term1-2005 2 Sound waves Why? � Are sound waves important? � Try solution x(p,t)=x 0 (p-ct) • Visually? Usually not � And x(p,t)=x 0 (p+ct) E � However, since speed of sound is a material � So speed of “sound” in rod is property, it can help us get to higher dimensions � � Speed of sound in terms of one spring (using � Courant-Friedrichs-Levy (CFL) condition: linear density m/L) is kL • Numerical methods only will work if information c = transmitted numerically at least as fast as in reality m (here: the speed of sound) � So in higher dimensions, just pick k so that c is • Usually the same as stability limit for good explicit constant methods [what are the eigenvalues here] • m is mass around spring [triangles, tets] • Implicit methods transmit information infinitely fast • Optional reading: van Gelder cs533d-term1-2005 3 cs533d-term1-2005 4
Potential energy Spring potential � Another way to look at the elastic spring forces: � For a single spring, through potential energy 2 � � x i � x j W ij = 1 � 1 � � � Recall for a system at position(s) x, potential 2 k ij L ij L ij � � energy W(x) gives the force F = � � W • Note we � re squaring the percent deformation (so this � x always increases as we move away from undeformed), and scaling by the strength of spring � For example, this gives conservation of total and by the length (amount of material it represents) dt + � W ( ) = v T M dv d dx energy K+E: 2 v T Mv + W 1 � To get the force on i, differentiate w.r.t. x i : � x dt dt � � x i � x j f ij = � � W ij x i � x j = v T Ma � F T v = � k ij � 1 � � � x i x i � x j � L ij � = v T F � v T F = 0 cs533d-term1-2005 5 cs533d-term1-2005 6 Directional derivatives 1D Continuum potential (regular calculus) � Add up the potential energies for each spring to � Pick a direction, or test vector, q approximate the total potential energy for the � Direction derivative along q is: elastic rod: 2 W ( x + � q ) � W ( x ) � � x i + 1 � x i � D q W ( x ) = lim W � � 1 ( p i + 1 � p i ) 1 2 E i + 12 � � � p i + 1 � p i � � � � 0 i � Take the limit as � p goes to zero: � Or alternatively 2 2 E ( p ) � x � � 1 D q W ( x ) = � g ( � ) = W ( x + � q ) g (0) where � W = � p � 1 1 dp � � � � � And the gradient � W/ � x is the vector s.t. 0 D q W ( x ) = � W � q � x i q � Now: how do we get forces out of this? The negative gradient of W w.r.t. x(p)? cs533d-term1-2005 7 cs533d-term1-2005 8
The variational derivative The variational derivative 2 � Now differentiate w.r.t. the scalar: � We want to take the “gradient” of a continuum 2 2 E � x � � p + � � q � 1 potential energy to get the force density: g ( � ) = d � � � p � 1 1 dp � � 2 f = � � W [ x ] 2 E � x � � 1 d � � � � W [ x ] = � p � 1 1 where dp � � 0 � x � � � � 2 � � x � � � x � � q � p + � 2 � q 1 2 0 2 E d � = � p � 1 + 2 � � p � 1 � But how do you differentiate w.r.t. a function 1 � � dp � � � � d � � � � � � p � � x(p)? 0 � Let � s first look at a directional derivative: look at � � � � x � � q � p + � � q 1 2 � the energy at x+ � q = � p � 1 E � dp � � � � � � p � � • q is a direction, or test function 0 2 2 E � x � � p + � � q � 1 � And evaluate at 0 to get the directional derivative � g ( � ) = W [ x + � q ] = � p � 1 E � x � � � q 1 1 dp � � � D q W [ x ] = g (0) = � � p � 1 � � � p dp � � � � 0 cs533d-term1-2005 9 0 cs533d-term1-2005 10 The variational derivative 3 Force density � The elastic force density is the negative of the � We want to make it look like an inner-product of variational derivative at a point: the “gradient” with q() • Use integration by parts: � � f = � � p E � x � � � p � 1 � � � � 1 � � � � � � � p E � x � � � + E � x � � � � 1 � D q W [ x ] = � � p � 1 � p � 1 � � � qdp � q � � � � � � � � � � � � � The acceleration of that point is force density 0 0 divided by mass density: � Ignoring boundary conditions for now, we see that the variational derivative of W at some � � � � � 2 x � p E � x � � t 2 = 1 � p � 1 interior point p is just: � � � � � � � � � � � � � � p E � x � � � p � 1 � � � � � Which is exactly what we got before! � � � � cs533d-term1-2005 11 cs533d-term1-2005 12
Discretized potential Multi-dimensional spring potential � Change the L scale in the 1D spring potential to � Now we have an alternative way to be the area/volume around the spring discretize the equation: • When we add up, we get an approximation of an • Approximate the potential energy integral with integral over the elastic object a discrete sum 2 � � x i � x j • Take the gradient to get forces W ij = 1 � 1 � � 2 E A ij L ij � � � This approach generalizes to all sorts of forces � Then get the spring force on i due to j: � � � Let � s do it for multi-dimensional springs x i � x j f ij = � � W ij x i � x j = � E A ij � 1 � � � x i x i � x j L ij � L ij � cs533d-term1-2005 13 cs533d-term1-2005 14 Bending Bending energy � Integrate mean curvature squared over the � For simulating cloth and thin shells, also surface: need forces to resist bending �� W bend = 2 B � 2 1 � Nontrivial to directly figure out such a force � on a triangle mesh � Let � s bypass the continuum formulation, and jump to a discrete approximation of this integral � Even harder to make sure it � s roughly � Split mesh up into regions around “hinges” mesh-independent (common edges between triangles) � So let � s attack the problem with a potential � W bend � 2 B e � e 2 A e 1 energy formulation e :edge • At each interior edge e, have bending stiffness B e , curvature estimate � e and area of region A e cs533d-term1-2005 15 cs533d-term1-2005 16
Edge curvature estimate Discrete bending energy � As a rough approximation then, using � Look at how the hinge is bent A e =|e|(h 1 +h 2 )/6: • Dihedral angle between the incident triangles 2 � � � e � � Think of fitting a cylinder of radius R parallel to W bend � 1 2 B e A e � � edge, mean curvature is 1/(2R) � � 6 A e e :edge • [side-view of triangle pair: angle between normals is theta, triangle altitudes are h 1 and h 2 ] 2 e � = � 2 1 2 B e � For small bend angles (limit as mesh is refined!) 36 A e � e :edge 1 2 R � � Treat all terms except theta as constant h 1 + h 2 (measured from initial mesh), differentiate that w.r.t. positions x cs533d-term1-2005 17 cs533d-term1-2005 18 Bending forces � See e.g. Bridson et al. “Simualtion of clothing…” for reasonable expressions for forces • Additional simplification: replace theta with a similarly-behaved trig function that can be directly computed • Warning: derivation is quite different! cs533d-term1-2005 19
Recommend
More recommend