Linear Subspace Exploration for Mechanical Systems

Project Members: Al Rahim Hossain, Stephanie Jung, Aleksa Milovanović, Reid Tang

Mentors: Otman Benchekroun, Ty Trusty

Introduction

Simulating deformable objects often requires optimizing over a large number of unknowns. For example, a 2D mesh with $n$ vertices has $2n$ positional variables: an $x$ and $y$ coordinate for each vertex. A 3D mesh similarly has $3n$ positional variables.

In a full order model, the deformed shape is represented by a vector $\mathbf{x}$ containing all vertex positions. It is found by minimizing the total energy of the object:

$\mathbf{x}^*=\underset{\mathbf{x}}{\text{argmin}}\left(E_{\mathrm{elastic}}(\mathbf{x})+E_{\mathrm{ext}}(\mathbf{x})\right)$

where $E_{\mathrm{elastic}}$ measures how much the object resists being deformed from its rest shape and $E_{\mathrm{ext}}$ represents external forces. The solution $\mathbf{x}^*$ is the state with the lowest total energy, where the elastic response and external forces are balanced. However, this solution may be computationally expensive as it treats every vertex position as an independently moving variable.

**Figure 1:** Comparison of the original 2D mesh (*left*) and its Neo-Hookean deformation (*right*) under gravity ( $g = -0.1$ ), with fixed vertices highlighted in red. This represents the exact solution $\mathbf{x}^*$ found using Newton’s method.

Reduced-Order Model

Reduced-order modelling reduces the size of this optimization by restricting the solution to a low-dimensional subspace. Instead of tracking every point independently, we choose a small set of meaningful deformation patterns. We then approximate the full configuration as $\mathbf{x} \approx \mathbf{x_0 + Bz}$ where $\mathbf{x_0}$ is the rest state of the object, the columns of matrix $\mathbf{B}$ represent the selected deformation patterns, and vector $\mathbf{z}$ contains their corresponding coefficients [1]. We can now solve for

$\mathbf{z}^* = \underset{\mathbf{z}}{\text{argmin }}E(\mathbf{x_0 + Bz})$

and the final full-space solution is given by $\mathbf{x}^* \approx \mathbf{x_0 + Bz}^*$ .

The effectiveness of a reduced-order model depends heavily on the choice of the basis $\mathbf{B}$ . Ideally, its columns should capture the important deformation patterns of the object using as few modes as possible. There are several ways to construct such a basis.

Proper Orthogonal Decomposition

Proper Orthogonal Decomposition constructs the basis from a dataset of observed or simulated deformations [2]. It finds the directions that best capture the variation present in the example shapes, making it a data-driven approach. It is found by solving

$\mathbf{B} = \underset{\mathbf{B}}{\text{argmin }} ||\mathbf{X – BB}^T\mathbf{X}||^2$ such that $\mathbf{B}^T\mathbf{B=I}$

Here, $\mathbf{X}$ is a matrix whose columns contain example displacement vectors from the rest configuration of the object. We use singular value decomposition, $\mathbf{X} = \mathbf{USV}^T$ , to find the deformation patterns that best represent these examples. The columns of $\mathbf{U}$ are ordered from most to least important so if we want a reduced space with 20 modes, for example, we simply take the first 20 columns $\mathbf{B} = \mathbf{U}[:, :20]$ .

**Figure 2:** *Top*: Full-order solutions with gravity $g \in [-0.002, -0.003, -0.004, -0.005, -0.006, -0.007]$ used to form the dataset. *Bottom*: Deformation under gravity $g = -0.1$ using POD. Fixed vertices highlighted in red.

A main drawback of POD is that it requires collecting a representative set of deformation examples beforehand. The resulting subspace can only represent deformation patterns present in, or similar to, this training data, so unseen motions may be captured poorly.

**Figure 3:** Deformation under gravity $g = -0.1$ using the POD basis from Figure 2 but with fixed vertices on the opposite ear. It fails since this configuration is not represented in the basis.

Linear Modal Analysis

Linear Modal Analysis builds the reduced basis from the object’s natural vibration modes around its rest shape [3]. The basis vectors are deformation modes obtained from the mesh’s stiffness and mass matrices, so the subspace is determined by the physics of the object rather than by example deformation data. These modes are found by solving the generalized eigenvalue problem

$\mathbf{H}\phi = \lambda\mathbf{M}\phi$ .

Here, $\mathbf{H}$ is the elastic stiffness of the object and $\mathbf{M}$ accounts for its mass. Solving the eigenvalue problem gives the natural deformation patterns $\phi$ , or modes, of the object. The lowest-frequency modes are usually the most important and are chosen as the columns of the reduced basis $\mathbf{B}$ .

**Figure 4:** *Top*: Animation of the first 10 deformation modes. *Bottom*: Deformation under gravity $g = -0.1$ using an LMA basis with 5 modes (left), 10 modes (center), and 20 modes (right).

However, LMA is based on a linearization around the rest shape, so it works best for small deformations. Large rotations or strongly nonlinear deformations may not be represented well by a basis of linear vibration modes.

Actuation-Aware Subspaces

POD can capture large deformations, but it requires representative simulation data. LMA avoids this data collection, but its basis is not aware of the parameters driving the simulation. Actuation-aware subspaces build the reduced basis directly from how the object deforms as the actuation parameters change. We write

$\mathbf{x(a)} = \underset{ \mathbf{x} } { \text{argmin } } E(\mathbf{x, a})$

where $\mathbf{a}$ may represent quantities such as spring rest lengths, gravity, muscle activation, or collider parameters. The goal is to approximate how $\mathbf{x(a)}$ changes as $\mathbf{a}$ changes.

We can first approximate this relationship using a first-order Taylor expansion around a reference actuation $\mathbf{a}_0$ :

$\mathbf{x}(\mathbf{a})\approx\mathbf{x}(\mathbf{a}_0)+\sum_i(a_i-a_{0,i})\frac{\partial\mathbf{x}}{\partial a_i}$ .

The derivatives $\frac{\partial \mathbf{x}}{\partial a_i}$ describe the deformation caused by changing each actuation parameter. For example, if the object is driven by springs, each parameter could control the rest length of one spring. However, this is still a linear approximation and works best for small deformations near $\mathbf{a}_0$ . Larger motions such as bending, rotation, and compression follow curved deformation paths that cannot be represented well using only fixed linear directions.

To capture these effects, we include second-order derivatives, $\frac{\partial^2 \mathbf{x}} {\partial a_i \partial a_j}$ , for the second-order Taylor expansion:

$\mathbf{x}(\mathbf{a}) \approx \mathbf{x}(\mathbf{a}_0) + \sum_i (a_i-a_{0,i}) \frac{\partial \mathbf{x}}{\partial a_i} + \frac{1}{2} \sum_i\sum_j (a_i-a_{0,i})(a_j-a_{0,j}) \frac{\partial^2 \mathbf{x}}{\partial a_i \partial a_j}$

If the first derivatives describe the response to actuation, the second derivatives describe how that response changes as the actuation changes. These terms help capture deformation patterns that appear together during larger motions, allowing the reduced space to follow nonlinear deformation paths more closely.

The reduced subspace can therefore be formed by

$\mathbf{B} = \left[ \frac{\partial \mathbf{x}}{\partial a_1} \ \cdots\ \frac{\partial \mathbf{x}}{\partial a_p} \ \frac{\partial^2 \mathbf{x}}{\partial a_1^2} \ \cdots\ \frac{\partial^2 \mathbf{x}}{\partial a_i \partial a_j} \right]$

such that this space is tailored to the parameters driving the simulation.

This idea is similar to modal derivatives [1], which describe how linear vibration modes change and interact. Here, we apply the same idea to actuation parameters, so the basis is tailored to the forces that will actually drive the simulation.

**Figure 5:** Deformation under spring actuation as the spring is made tighter using (from left to right): first-order approximation, second-order approximation, LMA, and an exact solution using Newton’s method.

**Figure 6:** Deformation error as the spring contracts for the first-order, second-order, and LMA reduced models compared with the exact solution.

Conclusion

Reduced-order models can make deformable simulations more efficient by representing motion with a small number of deformation patterns. In this project, we explored an actuation-aware basis that is built from how the equilibrium shape changes with the parameters driving the simulation. First-order derivatives capture the main response to actuation, while second-order derivatives help represent larger and more nonlinear deformations.

There are still several directions for future work. Hyper-reduction is needed to further reduce computational cost [4], and the current method relies on a PSD-projected Hessian for stable basis construction. The number of basis terms also grows with the number of actuation parameters. Future work could address these limitations and extend the method to more complex 3D scenes and different types of actuation.

References

[1] J. Barbič and D. L. James. Real-Time Subspace Integration for St. Venant-Kirchhoff Deformable Models. ACM Transactions on Graphics, 2005.

[2] L. Sirovich. Turbulence and the Dynamics of Coherent Structures. I. Coherent Structures. Quarterly of Applied Mathematics, 45(3), 561–571, 1987.

[3] A. Pentland and J. Williams. Good Vibrations: Modal Dynamics for Graphics and Animation. SIGGRAPH, 1989.

[4] C. Brandt, E. Eisemann, and K. Hildebrandt. Hyper-Reduced Projective Dynamics. ACM Transactions on Graphics, 2018.

Authors