# Physikalische Grundlagen für Computerspiele

Maximilian Stark (mail@dakror.de), WS2018

Stand: 17.02.2019

## Context

### Real-time physics

• Interactive -> Simulation
• Accuracy vs. Speed
• Stability
• Applications
• Training
• Education
• Entertainment

1. Geometric
2. Kinematic
3. Physical
4. Behavioral
5. Cognitive

### Animation techniques

• Conventional
• keyframe interpolation
• Full control
• Tedious
• No interaction
• Motion capture
• recording of markers on actor
• Real performance
• Bound to realism
• No real interaction
• Physics-based
• fully interactive
• Computationally expensive

### Goals

• Immersion
• Gameplay
• Exploration
• Asset creation reduction

## Mass spring systems

• Simulation of mesh structure by mass points with springs
• Hooke's law
• Spring deformation force proportional to deformation distance
• $\displaystyle F = -k(l-L)$, $\displaystyle k$ Stiffness, $\displaystyle l$ current length, $\displaystyle L$, original length
• Equations of motion (Newton 2nd law): $\displaystyle m_i*a_i = F_i^{\text{int}}(x_i(t)) + F_i^{\text{ext}}(t) - \underbrace{\gamma v_i(t)}_{\text{Damping}}$
• Implementation
• Euler step
• Generic force generators
• linked to single mass point
• evaluation & accumulation of force
• + Clean interface
• - Two spring instances per spring needed
• Simulation
• Analytic solution
• ODE: $\displaystyle x''=-\frac{k}{m}x$
• Harmonic oscillation with damping: $\displaystyle x(t) = e^{-\gamma t}d\cos(\sqrt{\frac{k}{m}}t)$
• - Only simplified setup
• Numerical integration
• From Taylor Series
• - Explosions
• Leapfrog
• 1 force evaluation, 2nd order accuracy
• Interleaved update of velocity and position
• $\displaystyle v(t+\frac{h}{2}) = v(t-\frac{h}{2}) + ha(t)$
$\displaystyle x(t+h) = x(t) + hv(t+\frac{h}{2})$
• Initialisation of velocity at $\displaystyle \frac{h}{2}$

### Numerical integration

• Reason
• $a(t) = v'(t) = x''(t)$
• Methods
• Explicit euler
• $\displaystyle y_{n+1} = y_n + hf(x_n,y_n)$
• Step in direction of current slope
• Heun
• $\displaystyle \widetilde{y} = y_0 + hf(x_0,y_0)$
$y_1 = y_0 + \frac{h}{2}(f(x_0,y_0) + f(x_1,\widetilde{y}))$
• Step with avg. tangent from current and from euler step → 2nd order
• "Predictor-corrector" method
• Midpoint
• $\displaystyle \widetilde{y} = y_0 + \frac{h}{2}f(x_0,y_0)$
$y_1 = y_0 + hf(x_0 + h/2,\widetilde{y})$
• derivative at half step
• Runge Kutta 2nd order
• Runge Kutta
1. First estimate of slope: $\displaystyle k_1 = f(x_0,y_0)$
2. Prediction of tangent at midpoint: $\displaystyle k_2 = f(x_0 + \frac{h}{2}, y_0 +\frac{h}{2} * k_1)$
3. Correction of prediction: $\displaystyle k_3 = f(x_0 + \frac{h}{2}, y_0 + \frac{h}{2} * k_2)$
4. Prediction of slope at full step $\displaystyle k_4 = f(x_0 + h, y_0 + h*k_3)$
5. Weighted step taken: $\displaystyle y_1 = y_0 + \frac{h}{6} (k_1 + 2k_2+2k_3+k_4)$
• Velocity verlet
• Higher order interpolation instead of multiple evaultions of derivative
• $\displaystyle y_{n+1} = y_n + hy_n' + \frac{h^2}{2}y_n''+\mathcal{O}(h^3)$
• $\displaystyle y_{n+1} = y_n' + h \frac{y_n'' + y_{n+1}''}{2} + \mathcal{O}(h^3)$
• data vs. calculation tradeoff ($y_n''$ stored)
• Implicit integration
• Approximation of derivative at next time step instead of current
• $\displaystyle y_{n+1} = y_n + hf(x_{n+1}, y_{n+1})$
• Solution of (linear) $\displaystyle Ay_{n+1}=b$
• + Unconditionally stable
• - High complexity
• - Uncontrollability when non-linear
• - Boring aesthetics
• Accuracy
• n step 1st-order much worse than 1 step n-order
• Speed vs accuracy → Visual validity
• Error bounded by n-th derivative in interval: $\displaystyle 0 \leq e \leq \frac{h^n}{n}y^n$

## Rigid Bodies

• Center of mass $\displaystyle x_{cm} = \frac{\sum_i m_i x_i}{\sum_i m_i}$
• Angular velocity $\displaystyle w$
• Angular momentum $\displaystyle L = \sum_i (x_i \times m_i v_i)$
• Torque $\displaystyle q = \frac{\mathrm{d}L}{\mathrm{d}t} = \sum_i (x_i \times f_i)$: Rotational force
• Inertia tensor $\displaystyle i$: Resistance to rotation

### 2D

• Inertia tensor scalar: $\displaystyle i = \sum_n (m_n x_n \cdot x_n)$
• Change of angular velocity: $\displaystyle w(t+h) = w(t) + h\frac{q}{i}$
• Change of rotation: $\displaystyle r(t+h) = r(t) + hw(t)$
• Simulation
1. External forces
2. Euler Step
1. Position
2. Velocity
3. Rotation
4. Angular velocity
3. World position
• Collision detection
• Relative velocity: $\displaystyle v_{rel} = n \cdot (v_A - v_B)$
• $\displaystyle v_{rel} < 0$: Collision
• $\displaystyle v_{rel} > 0$: Separation
• $\displaystyle v_{rel} = 0$: Sliding
• Collision response
• using forces
• unpredictable
• full rigidity → infinite force
• Impulses
• Velocity change without time
• $\displaystyle J = m\Delta v$
• linear motion: $\displaystyle J_{lin} = \frac{-(1+c)v_{rel}\cdot n}{\frac{1}{M_a}+\frac{1}{M_b}}$
• full motion: $\displaystyle J = \frac{-(1+c)v_{rel}\cdot n}{\frac{1}{M_a}+\frac{1}{M_b} + \frac{(x_a \times n)^2}{I_a} + \frac{(x_b \times n)^2}{I_b}}$
• $\displaystyle v_{rel} < 0$ → $\displaystyle J >= 0$ scalar
• $\displaystyle v_a' = v_a + \frac{Jn}{M_a}$
• $\displaystyle v_b' = v_b - \frac{Jn}{M_b}$
• $\displaystyle w_a' = w_a + \frac{(x_a \times Jn)}{I_a}$
• $\displaystyle w_b' = w_b + \frac{(x_b \times Jn)}{I_b}$
• Bounciness: coefficient of restitution
• $\displaystyle c=1$: fully elastic
• $\displaystyle c=0$: fully plastic
• Loss of kinetic energy for $\displaystyle c<1$ → deformation of body

### Rotations

• Fixed Angles
• Reverse order: $\displaystyle R_{x,y,z} = R_z R_y R_x$
• Euler Angles
• Axes fixed to object
• Yaw (z), Pitch (x), Roll (y)
• Quaternions
• $\displaystyle q = (\cos(\frac{\phi}{2}), n\sin(\frac{\phi}{2}))$, angle $\displaystyle \phi$, axis $\displaystyle n$
• Rotation of point $\displaystyle p$: $\displaystyle q(0,p)\widetilde{q}=p'$ with $\displaystyle \widetilde{q} = (s, -v)$
• $\displaystyle q_1 \cdot q_2 = (s_1, v_1) \cdot (s_2, v_2) = (s_1 s_2 - v_1 \cdot v_2, ~ s_1 v_2 + s_2 * v_1 + v_1 \times v_2)$

### 3D

• Inertia Tensor -$3 \times 3$ matrix, Computation using mass-weighted co-variance matrix of body coordinate positions
$I = \mathrm{Id} ~ \mathrm{trace}(C)-C; ~ C = \sum_n m_n x_n x_n^\intercal; ~ C_{j,k} = \sum_n m_n, x_{n,j} x_{n,k}$
• Dependant on current orientation
• Generally working with inverse
• Computation from initial values: $\displaystyle I_\text{current} = \mathrm{Rot}_r ~I_0~ \mathrm{Rot}_r^{-1} = \mathrm{Rot}_r ~I_0~ \mathrm{Rot}_r^\intercal$
• Angular motion
• velocity not constant, momentum constant
• velocity dependant on rotational axis
• $\displaystyle L=Iw$
• $\displaystyle w(t+h) = I^{-1} ~ L(t+h)$
• Orientation
• Derivative of quaternion: $\displaystyle \frac{\mathrm{d}r}{\mathrm{d}t} = \frac{1}{2}\binom{0}{w}$
• Integration using $\displaystyle r' = r + \frac{h}{2} \binom{0}{w} r$
• Rotation matrix integration using shear matrix
• Simulation
1. External forces
2. Euler Step
1. Position
2. Velocity
3. Rotation
4. Angular velocity
5. Inverse inertia tensor
6. Angular momentum
3. World position
• Inapplicability of leapfrog: Forces dependant on position
• Impulse
• $\displaystyle J = \frac{-(1+c)v_{rel}\cdot n}{\frac{1}{M_a}+\frac{1}{M_b} + \left[ (I_a^{-1}(x_a \times n)) \times x_a + (I_b^{-1}(x_b \times n)) \times x_b \right] \cdot n}$
• $\displaystyle v_a' = v_a + \frac{Jn}{M_a}$
• $\displaystyle v_b' = v_b - \frac{Jn}{M_b}$
• $\displaystyle L_a' = L_a + (x_a \times Jn)$
• $\displaystyle L_b' = L_b - (x_b \times Jn)$

## Fluid Simulation

• Velocity $\displaystyle u$
• Pressure $\displaystyle p$
• Density $\displaystyle \rho$
• Conservation of momentum: $\displaystyle \text{advection} = \text{pressure} + \text{viscosity} + g$
• Conservation of mass
• Divergence: Existance of sink / source

### Representations

• Eulerian: grid based property field
• + Precision
• + Ease of calculation
• - Problems with continuous motion
• Conservation of mass: mass per cell, sum constant
• Lagrangian: particle like

### Smoothed Particle Hydrodynamics (SPH)

• Solution of NS equations using particles
• Euler integration of position & velocity
• Goal: Calculation of forces for NS-equation approximation
• Smooth particles
• Mass as continuous field
• Smoothing kernel $\displaystyle W$ for interpolation
• Attribute evaluation: $\displaystyle a_W(x) = \sum_i a_i \frac{m}{\rho_i}W(x-x_i)$
• Gradient evaluation: $\displaystyle \nabla a_W(x) = \sum_i a_i \frac{m}{\rho_i} \nabla W(x-x_i)$
• Forces
• Pressure
• Conversion of density using equation of state
$p_i = k((\frac{\rho_i}{\rho_0})^7 - 1); ~ k$ Stiffness
• $\displaystyle f_\text{pressure} (x_j) = - \sum_i \frac{p_i + p_j}{2} \frac{m}{\rho_i} \nabla W$
• Viscosity
• $v$ viscosity, internal friction
• $\displaystyle f_\text{visc}(x_j) = - \sum_i (v_i - v_j) \frac{m}{\rho_i} \nabla^2 W$
• Color field
• $\displaystyle c(x_j) = - \sum_i \frac{m}{\rho_i} W$
• Surface tension
• $\sigma$
• $c$
• $\displaystyle f_\text{st} = -\sigma \nabla^2 c \frac{n}{|n|}$
• Simulation
1. Computation of grid / kd-tree
2. Density & color-field
3. Forces
1. Pressure
2. Viscosity
3. Surface tension
4. Integration of position & velocity (Euler)
• Issues
• Collision with obstacles: Clumping
• No flow orthogonal to walls
• High pressure forces from clusters
• Negative pressure: single pressure smaller than of rest
• Explicit integration: stability problems
• (Too) many interacting parameters & systems
• Improvement
• Position based method
• Velocity based on position update
• Predictions
• No iterative force computation

#### Smoothing kernels

• radius $\displaystyle d$
• Functions
• Cubic splines, distance to center $\displaystyle q$
• Polynomial functions
• Key properties
• Fast evaluation
• Computation of derivatives
• Improvements
• Different kernels for different forces
• Kernel with non-negative 2nd derivative
• Spiky kernel: high values near origin

## Collision detection & response

#### Culling

• Sweep and Prune
• SAT with AABBs
• $\mathcal{O}(n\log n + k)$: sorting ($\mathcal{O}(n\log n)$) + reporting intersection: $\mathcal{O}(k)$
• Sorting with temporal coherence (insertion sort): $\mathcal{O}(n)$
• Cullide
• Visibility based
• No collision if fully in front of them (z-buffer pixel count)
• $\mathcal{O}(n)$
• Profits from hardware acceleration

#### Bounding Volume Hierarchies

• Bounding Volumes (asc quality, desc complexity)
• Spheres
• AABBs
• OBBs
• axes: Eigenvectors of covariance matrix $C$
• $\displaystyle \mu = \frac{1}{n}\sum_i x_i; ~ C_{j,k} = \frac{1}{n} \sum_i (x_{i,j}-\mu_j)(x_{i,k} - \mu_k)$
• SAT using 15 axes: 3x face normals of A, 3x face normals of B, 3*3 edge pairs
• k-discrete oriented polytopes (k-DOP)
• k fixed directions
• k max values
• SAT on k axes
• Convex hulls
• Hierarchies
• Subdivision of volumes to generate hierarchy
• Bottom-up: from triangles to hierarchy by grouping
• Top-down: recursive splitting
• Insertion
• Graph building by proximity grouping
• - Quality dependant on insertion order
• Collision tests
• SAT
• Voronoi marching
• GJK algorithm
• Deformable objects
• Sphere tree
• Update of leaves only for deformed primitives
• Information propagation to root
• + Efficient
• - Suboptimal BV

#### Spatial partitioning

• Data structures
• Uniform grid
• kd-tree
• BSP-tree
• Binary space partitioning tree
• Generalization of kd-tree
• Creation
• Selection of plane with arbitrary normal
• Sorting of geometry into the two spaces
• Splitting of polygons intersection selected plane
• Recursive division
• - time consuming
• Useful for efficient rendering
• Spatial hashing for self-collisions
• Heavily twisted deformable objects
• Volumetric meshes
• tetrahedal mesh
• Hash function $H(i,j,k) = (ip_1 \oplus jp_2 \oplus kp_3) ~ \mathrm{mod} ~ n$
• $i,j,k$: cell coordinates
• $p_1,p_2,p_3$: large primes
• $n$: hash table size
• Cell based partition & hashing
• Implicit uniform grid
1. Vertex hashing: {cell -> list of vertices in cell}
2. Tetraheda hashing
• tetraheda hashed to all cells touching bounding box
3. Vertices and tetraheda in same cells tested for intersection
• Vertex-Tetrahedon collision test
• Barycentric coordinate based
• All points $P$ in tetrahedon: $P = (b_1, b_2, b_3); ~ b_1+b_2+b_3 = 1; ~ b_i = a_i / (a_1+a_2+a_3)$
• Intersection if: $\mathrm{coord}(A) = \mathrm{vol}(PBCD) / \mathrm{vol}(ABCD) \in (0,1) ~ \forall$ vertices
• Oriented faces
• $n$ face normal
• $n \cdot (P-A) > 0 ~ \forall$ faces
• Signed-distance functions
• Sampling of distance to object with dense grid
• Indication of inside (neg.) / outside (pos.) by sign of value
• part of level set framework (Niveaumengen)
• Embedding of surface into higher dimensional representation
• flat surface -> volumetric function
• Surface given by level set $\Gamma = \{ x | \phi (x)=0\}; ~ \phi: \mathbb{R}^3\mapsto \mathbb{R}$
• Implicit modification of surface
• Calculatation of Closest Surface Point
• Gradient always pointing away from object (Normalization)
• Data Acquisition
• Precomputation
• Triangle rasterization
• Ray intersections
• Data from previous stages
• Pros/Cons
• + Fast checking
• + Single precomputation
• - Lots of memory
• Improvements using adaptive grid structure

### Narrow Phase

#### Contact & Friction

• Stationary objects
• Balance of forces
• Gravity canceled by reaction force from contact
• → Resting contacts
• Multiple impulses
• Change of total velocity by every impulse
• Final relative velocity must be separating:$\Delta v_i + \sum_k J_k \geq 0$
• Solution in single system of equations: No unique solution
• → Constraint: Impulses always positive and minimal
• Linear system matrix as representation of collision influence graph
• Micro Collisions
• Resting contacts as small collisions
• Removal of gravity acceleration from last time step
• Detection of tiny velocities, decrease of restitution coefficient
• → Cancels out bounce back
• Impulse without gravity
• Friction
• static: material parameter, solid ground
• dynamic: two moving, touching objects
• Impulse with friction: apply impulse & friction (along tangent)

#### Position correction

• Approaches
• Teleport to surface of partner
• Backtrace velocity & orientation
• Account for mass
• Minimal movement
• Based on penetration depth instead of $v_{\text{rel}}$
• Orientation change: possibility of new collisions
• minimize change
• e.g proportional to rel. collision point on object

## Full RB simulation

1. Integration
2. Collision detection
• Obtaining list of colliding object pairs
• Usage of accelerated data structure
• Minimal face-corner collisions
• List building over time
1. Timestep
• Gather single deepest collections
2. Timestep
• Collision verification
• Check of $v_{\text{rel}}$
3. Penetration resolution
• Sorted list, descending
• Impulse based PosRot change
• In theory: re-check for collisions
• Usually finite time / attempts for resolution
• Iterative Convergence
• Priority based
4. Velocity resolution
• Sorted list, ascending $v_{\text{rel}}$
• Recomputation of velocity for all objects colliding with either partner
• Collision tolerance
• Ignoring of small collisions
• Minimal penetration depth & velocity
• Only visual validity

### Sleeping

• Performance boost
• Sleep state basis
• kinetic energy
• velocity
• recent history (EMA: exponential moving average)
• Still participation in CD with 0 velocity
• Waking up
• Collision with active object
• Hysteresis (threshold range)

## Complex cloth simulation

• Mass spring system with springs for stretching & bending
• Thickness $c$
• Difficulty of tangling prevention
• Collision detection
• AABB Hierarchy
• Bounding box for movement
• Including old & new position
• Thickness
• Safety margin
• Visual validity
• Fully plastic collisions
• Position correction: relaxation (factor 0.25)
• Penetration depth = thickness → Collision ignored, friction-esque
• Heun integration
• Smoothed rendering
• Visual collision checking based on vertex displacement "velocity"
• Temporal incoherence