Title: Position Based Dynamics

Authors: Matthias Müller, Bruno Heidelberger, Marcus Hennix, John Ratcliff

Published in: VRIPHYS 07: Fourth Workshop on Virtual Reality Interactions and Physical Simulations (2006)

• 해당 포스팅은 AI를 통해 번역한 내용입니다. 일부 내용이나 해석에 오역 또는 의역이 있을 수도 있습니다.
• 개인적으로 공부하며 핵심 내용을 정리하고 해설한 포스팅입니다.
• 본 저작물의 원문 저작권은 저자 및 출판사에 있습니다. 문제가 될 시 삭제될 수 있으며, 원문은 영리적 목적으로 복사 혹은 배포가 금지되어 있음을 알립니다.

1. Introduction

컴퓨터 그래픽스 분야의 물리 기반 애니메이션 연구는 강체(rigid bodies), 변형 가능한 물체(deformable objects), 또는 유체 흐름(fluid flow)과 같은 물리적 현상을 시뮬레이션하기 위한 새로운 방법들을 찾는 데 주력하고 있습니다. 정확성에 주안점을 두는 계산 과학(computational sciences)과는 대조적으로, 이 분야의 주요 쟁점은 결과가 시각적으로 그럴듯해 보이면서도 안정성(stability), 견고함(robustness), 그리고 속도를 확보하는 것입니다.
따라서 계산 과학의 기존 방법들을 그대로 채택하여 사용할 수는 없습니다. 사실, 컴퓨터 그래픽스에서 물리 기반 시뮬레이션을 연구하는 주요 명분은 이 분야의 특수한 요구 사항에 맞춘 전문화된 방법을 고안해 내는 것입니다. 우리가 제시하는 방법 또한 이러한 범주에 속합니다.

동적 객체를 시뮬레이션하는 전통적인 접근 방식은 힘(force)을 이용하는 것이었습니다. 각 타임 스텝(time step)의 시작 부분에서 내부 및 외부의 힘들이 축적됩니다. 내부 힘의 예로는 변형 가능한 물체의 탄성력이나 유체의 점성 및 압력 등이 있으며, 중력과 충돌력은 외부 힘의 대표적인 예시입니다.
뉴턴의 운동 제2법칙은 질량을 매개로 힘과 가속도를 연결합니다. 따라서 밀도나 정점의 집중 질량(lumped masses)을 이용하여, 힘은 가속도로 변환됩니다. 그다음 어떠한 시간 적분법(time integration scheme)을 사용하더라도, 가속도로부터 속도를 먼저 계산하고 그 속도로부터 위치를 계산할 수 있습니다. 일부 접근 방식은 애니메이션을 제어하기 위해 힘 대신 충격량(impulses)을 사용하기도 하는데, 충격량은 속도를 직접적으로 변화시키기 때문에 적분 단계 하나를 생략할 수 있습니다.

컴퓨터 그래픽스, 특히 컴퓨터 게임에서는 물체의 위치나 메쉬(mesh)의 정점(vertex)을 직접 제어해야 하는 경우가 많습니다. 사용자는 특정 정점을 운동학적(kinematic) 물체에 부착하거나, 정점이 항상 충돌하는 물체의 바깥에 머물도록 확실히 조치하고 싶어 할 수 있습니다. 우리가 여기서 제안하는 방법은 위치를 직접 조작하기 때문에 이러한 제어를 쉽게 만들어 줍니다.
또한, 위치 기반 접근 방식을 사용하면 적분 과정을 직접 제어할 수 있어, 명시적 적분(explicit integration) 방식에서 발생하는 오버슈팅(overshooting)이나 에너지 증가 문제를 피할 수 있습니다. 요약하자면, 위치 기반 역학의 주요 특징과 장점은 다음과 같습니다.

• 위치 기반 시뮬레이션은 명시적 적분을 제어할 수 있게 해주며, 전형적인 불안정성 문제를 제거합니다.
• 시뮬레이션 도중에 정점의 위치나 물체의 일부를 직접 조작할 수 있습니다.
• 우리가 제안하는 공식화 방식은 위치 기반 환경에서 일반적인 제약 조건(constraints)을 처리할 수 있게 합니다.
• 명시적 위치 기반 솔버(solver)는 이해하기 쉽고 구현하기 용이합니다.

---

2. Related Work

최신 기술 보고서 [NMK∗05]는 질량-용수철 시스템, 유한 요소법(FEM), 또는 유한 차분법 등 변형 가능한 물체를 시뮬레이션하기 위해 컴퓨터 그래픽스 분야에서 사용되는 여러 방법들에 대해 훌륭한 개요를 제공하고 있습니다. 하지만 [MHTG05]의 인용 사례를 제외하면, 위치 기반 역학(position based dynamics)이라는 개념은 이 조사 보고서에 직접적으로 등장하지 않습니다. 그럼에도 불구하고, 위치 기반 접근 방식의 일부 요소들은 그 명칭을 명시하거나 완전한 프레임워크로 정의하지 않았을 뿐, 이미 다양한 논문에서 부분적으로 사용되어 왔습니다.

Jakobsen [Jak01]은 위치 기반 접근 방식을 바탕으로 자신의 Fysix 엔진을 구축했습니다. 그의 핵심 아이디어는 베를레 적분기(Verlet integrator)를 사용하여 위치를 직접 조작하는 것이었습니다. 이 방식에서는 속도가 현재 위치와 이전 위치 값 속에 암시적으로 포함되어 있으므로, 위치를 조작하면 속도 또한 자동으로 업데이트됩니다. 그는 주로 거리 제약 조건에 집중했으나, 더 일반적인 제약 조건들을 처리하는 방법에 대해서는 모호한 실마리만을 남겼습니다.
반면, 본 논문에서 우리는 일반적인 제약 조건들을 모두 처리할 수 있는 완전히 일반화된 접근 방식을 제시합니다. 또한 위치 투영(position projection) 과정에서 발생하는 선운동량 및 각운동량 보존이라는 중요한 이슈를 깊이 있게 다룹니다. 특히 이전 위치를 저장하는 방식 대신 명시적인 속도 값을 사용함으로써, 감쇠(damping)나 마찰(friction) 시뮬레이션을 훨씬 용이하게 구현할 수 있도록 했습니다.

Desbrun [DSB99]과 Provot [Pro95]는 질량-용수철 시스템에서 용수철이 과도하게 늘어나는 현상을 방지하기 위해 제약 조건 투영(constraint projection) 기법을 사용합니다. 하지만 이는 본 논문에서 제안하는 완전한 위치 기반 접근 방식과는 차이가 있습니다. 이들의 방식에서 투영은 지나치게 늘어난 일부 용수철을 다듬기 위한 후처리(polishing) 공정으로만 활용될 뿐, 시뮬레이션의 근간이 되는 핵심 방법론으로 채택된 것은 아닙니다.

Bridson 등은 천 시뮬레이션을 위해 전통적인 힘 기반 접근 방식을 사용하면서, 충돌 해결을 위한 충격량이 안정적인 범위 내에 머물도록 위치 기반의 기하학적 충돌 해결 알고리즘을 결합하여 사용합니다. 이는 Volino 등이 제안한 운동학적 충돌 수정 단계에서도 마찬가지로 나타나는 특징입니다.

Clavet 등 [CBP05]은 점탄성 유체(viscoelastic fluids)를 시뮬레이션하기 위해 위치 기반 접근 방식을 도입했습니다. 하지만 그들의 방식은 위치 투영 과정의 여러 단계에서 타임 스텝(시간 간격) 변수가 관여하고 있어, 진정한 의미의 완전한 위치 기반 방식이라고 하기는 어렵습니다. 따라서 이 적분법은 일반적인 명시적 적분(explicit integration)과 마찬가지로 조건부 안정성(conditionally stable)만을 가집니다.

Müller 등 [MHTG05]은 객체의 정지 상태(rest state)를 현재 상태(current state)에 매칭하여 찾아낸 특정 목표 위치(goal positions)로 점들을 이동시킴으로써 변형 가능한 물체(deformable objects)를 시뮬레이션합니다. 그들의 적분 방법(integration method)은 우리가 여기서 제안하는 방식과 가장 유사한 형태를 띠고 있습니다. 다만, 그들은 하나의 특수한 전역 제약 조건(global constraint)만을 처리하기 때문에, 본 논문에서 제안하는 것과 같은 별도의 위치 솔버(position solver)를 필요로 하지는 않습니다.

Fedor [Fed05]는 게임 내 캐릭터(characters)를 시뮬레이션하기 위해 Jakobsen의 접근 방식을 사용합니다. 그의 방법은 인간 캐릭터를 시뮬레이션하는 특정 문제에 최적화되어 있습니다. 그는 여러 개의 골격 표현(skeletal representations)을 활용하며, 투영(projections) 기법을 통해 이들 사이의 동기화(sync)를 유지합니다.

Faure [Fau98]는 속도(velocities)보다는 위치(positions)를 직접 수정하는 방식의 베를레 적분법(Verlet integration scheme)을 사용합니다. 하지만 우리가 비선형 제약 함수(non linear constraint functions)를 직접 다루는 것과 대조적으로, 그는 제약 조건(constraints)을 선형화(linearizing)하여 새로운 위치를 계산한다는 차이가 있습니다.

우리는 [BW98] 및 [THMG04]의 연구들처럼 제약 함수(constraint function)를 통해 일반적인 제약 조건들을 정의합니다. 이때 제약 함수 에너지의 도함수(derivative)를 구해 힘(forces)을 계산하는 기존 방식 대신, 평형 상태(equilibrium configuration)를 직접적으로 구해내고 위치를 투영(project positions)하는 방식을 취합니다. 이러한 방법을 통해 우리는 [GHDS03]과 [BMF03]에서 제안된 것과 유사하면서도 점 기반 접근 방식(point based approach)에 최적화된 천의 굽힘 항(bending term)을 도출해 냈습니다.

섹션 4에서는 이러한 위치 기반 역학(position based dynamics) 접근법을 활용한 천 시뮬레이션(simulation of cloth) 사례를 다룹니다. 천 시뮬레이션은 최근 컴퓨터 그래픽스 분야에서 매우 활발히 연구되어 온 영역으로, 개별 논문들을 일일이 나열하는 대신 종합적인 서베이 논문인 [NMK∗05]를 참고 문헌으로 제시합니다.

---

3. Position Based Simulation

이 섹션에서는 일반적인 위치 기반 접근 방식(position based approach)을 공식화하여 설명합니다. 이어지는 섹션과 결과 섹션에서는 본 방법론의 구체적인 응용 사례로서 천 시뮬레이션(cloth simulation)을 살펴볼 것입니다. 본 논문에서는 3차원 세계(three dimensional world)를 배경으로 고려하고 있으나, 제안하는 접근 방식은 2차원(two dimensions) 환경에서도 동일하게 훌륭하게 작동합니다.

3-1. Algorithm Overview

우리는 $N$개의 정점(vertices)과 $M$개의 제약 조건(constraints)의 집합으로 동적 객체(dynamic object)를 표현합니다.
정점 $i \in [1, \dots, N]$는 질량 $m_i$, 위치 $\mathbf{x}_i$, 그리고 속도 $\mathbf{v}_i$를 가집니다.

제약 조건 $j \in [1, \dots, M]$는 다음과 같이 구성됩니다.

• 카디널리티(cardinality, 연관된 정점의 수) $n_j$,
• 함수(function) $C_j : \mathbb{R}^{3n_j} \to \mathbb{R}$,
• 인덱스 집합(set of indices) $\{i_1, \dots, i_{n_j}\}, i_k \in [1, \dots, N]$,
• 강성 파라미터(stiffness parameter) $k_j \in [0 \dots 1]$
• 등식(equality) 또는 부등식(inequality) 중 하나의 타입.

등식(equality) 타입의 제약 조건 $j$는 $C_j(\mathbf{x}_{i_1}, \dots, \mathbf{x}_{i_{n_j}}) = 0$일 때 만족됩니다. 부등식(inequality) 타입일 경우, $C_j(\mathbf{x}_{i_1}, \dots, \mathbf{x}_{i_{n_j}}) \ge 0$이면 만족됩니다.
강성 파라미터(stiffness parameter) $k_j$는 0에서 1 사이의 범위에서 제약 조건의 강도를 정의합니다.
이 데이터와 타임 스텝(time step) $\Delta t$를 기반으로, 동적 객체는 다음과 같이 시뮬레이션됩니다.

(1) forall vertices $i$
(2) $\quad$ initialize $\mathbf{x}_i = \mathbf{x}_i^0, \mathbf{v}_i = \mathbf{v}_i^0, w_i = 1/m_i$
(3) endfor
(4) loop
(5) $\quad$ forall vertices $i$ do $\mathbf{v}_i \leftarrow \mathbf{v}_i + \Delta t w_i \mathbf{f}_{ext}(\mathbf{x}_i)$
(6) $\quad$ dampVelocities($\mathbf{v}_1, \dots, \mathbf{v}_N$)
(7) $\quad$ forall vertices $i$ do $\mathbf{p}_i \leftarrow \mathbf{x}_i + \Delta t \mathbf{v}_i$
(8) $\quad$ forall vertices $i$ do generateCollisionConstraints($\mathbf{x}_i \to \mathbf{p}_i$)
(9) $\quad$ loop $solverIterations$ times
(10) $\quad \quad$ projectConstraints($C_1, \dots, C_{M+M_{coll}}, \mathbf{p}_1, \dots, \mathbf{p}_N$)
(11) $\quad$ endloop
(12) $\quad$ forall vertices $i$
(13) $\quad \quad \mathbf{v}_i \leftarrow (\mathbf{p}_i - \mathbf{x}_i) / \Delta t$
(14) $\quad \quad \mathbf{x}_i \leftarrow \mathbf{p}_i$
(15) $\quad$ endfor
(16) $\quad$ velocityUpdate($\mathbf{v}_1, \dots, \mathbf{v}_N$)
(17) endloop

(1)-(3)번 라인은 시뮬레이션에 필요한 상태 변수(state variables)들을 초기화하는 단계입니다. 위치 기반 역학(position based dynamics)의 핵심 아이디어는 (7)번, (9)-(11)번, 그리고 (13)-(14)번 라인에 잘 나타나 있습니다.
우선 (7)번 라인에서는 명시적 오일러 적분 단계(explicit Euler integration step)를 통해 정점들이 이동할 새로운 위치의 예측치(estimates) $\mathbf{p}_i$를 계산합니다.
이후 (9)-(11)번 라인의 반복적 솔버(iterative solver)가 이 위치 예측치들이 제약 조건(constraints)을 만족하도록 조정합니다.
이 과정은 각 제약 조건을 가우스-자이델 방식(Gauss-Seidel type fashion)으로 반복하여 투영(project)함으로써 이루어집니다.
(13)번과 (14)번 단계에 이르면, 정점들의 실제 위치를 최적화된 예측치로 이동시키고 그에 따라 속도를 업데이트합니다. 이러한 방식은 현재 위치와 이전 위치의 차이를 통해 속도를 암시적으로 저장하는 베를레 적분 단계(Verlet integration step) 및 위치 수정 방식[Jak01]과 논리적으로 정확히 일치합니다.
그럼에도 불구하고, 본 논문처럼 속도(velocities) 데이터를 명시적으로 사용하는 방식은 물리량을 조작하는 데 있어 더욱 직관적인 통제 수단을 제공한다는 장점이 있습니다.

시뮬레이션 루프 내에서 속도(velocities)의 조작은 (5)번, (6)번, 그리고 (16)번 라인에서 이루어집니다. 우선 (5)번 라인은 위치 제약 조건(positional constraints)으로 직접 변환하기 어려운 외부 힘(external forces)이 있을 때 이를 시스템에 반영하는 역할을 합니다. 본 논문에서 우리는 이 단계를 중력을 추가하는 용도로만 사용하며, 이때 식은 $\mathbf{v}_i \leftarrow \mathbf{v}_i + \Delta t \mathbf{g}$($\mathbf{g}$는 중력 가속도)로 표현됩니다.
(6)번 라인에서는 시뮬레이션의 안정성을 위해 필요한 경우 속도 감쇠(damping)를 수행합니다. 특히 섹션 3.5에서는 물체의 전체적인 이동이나 회전과 같은 강체 모드(rigid body modes)를 방해하지 않으면서 전역 감쇠(global damping)를 적용하는 기법을 상세히 설명합니다. 마지막으로 (16)번 라인에서는 충돌이 발생한 정점들의 속도를 마찰(friction) 및 반발 계수(restitution coefficients)에 근거하여 최종적으로 수정하며 한 스텝의 계산을 마무리합니다.

시뮬레이션에서 제공되는 기본 제약 조건(constraints) $C_1, \dots, C_M$은 전체 과정 동안 고정(fixed)된 상태를 유지합니다. 하지만 이와 별개로, 알고리즘의 (8)번 라인에서는 매 타임 스텝(time step)마다 상황에 따라 변화하는 $M_{coll}$개의 충돌 제약 조건(collision constraints)을 실시간으로 생성합니다. 최종적으로 (10)번 라인의 투영 단계(projection step)에서는 이 고정된 제약 조건들과 가변적인 충돌 제약 조건들 모두를 통합하여 고려하고 처리합니다.

이 방식은 무조건적으로 안정적(unconditionally stable)인 특성을 가집니다. 이는 (13)번과 (14)번의 적분 단계가 전통적인 명시적 기법(explicit schemes)들처럼 미래의 상태를 맹목적으로 추정(extrapolate)하여 계산하지 않기 때문입니다. 대신, 제약 조건 솔버(constraint solver)를 통해 도출된 물리적으로 유효한 구성(configuration)인 $\mathbf{p}_i$로 정점들을 직접 이동시키는 방식을 취합니다.
물론 불안정성이 발생할 가능성이 아예 없는 것은 아니나, 그 원인은 오직 유효한 위치를 찾기 위해 뉴턴-랩슨 방법(Newton-Raphson method)을 사용하는 솔버 그 자체에 국한됩니다. 중요한 점은, 이 솔버의 안정성이 타임 스텝(time step)의 크기에 좌우되는 것이 아니라 제약 함수(constraint functions)가 어떤 형태를 띠고 있느냐에 따라 결정된다는 사실입니다.

이 적분법(integration)은 명시적(explicit) 또는 암시적(implicit) 기법(schemes)의 범주 중 어느 한쪽으로 명확하게 분류되지 않는 특성을 가집니다. 만약 타임 스텝(time step)당 솔버 반복(solver iteration)을 단 한 번만 수행한다면, 이는 명시적 기법(explicit scheme)과 상당히 유사한 모습을 보입니다. 그러나 반복 횟수를 늘려감에 따라 제약 조건이 있는 시스템(constrained system)을 원하는 만큼 임의로 뻣뻣하게(stiff) 만들 수 있으며, 이때 알고리즘은 암시적 기법(implicit scheme)에 더 가깝게 거동하게 됩니다. 이러한 반복 횟수의 증가는 시뮬레이션의 계산 병목(bottleneck) 지점을 충돌 감지(collision detection) 단계에서 솔버(solver) 단계로 전환시킵니다.

💡 왜 '위치'를 먼저 계산할까?
기존의 힘 기반(Force-based) 시뮬레이션은 힘 → 가속도 → 속도 → 위치 순으로 계산합니다. 하지만 이 방식은 힘이 조금만 커져도 에너지가 폭발하여 시스템이 깨지기 쉽습니다.
PBD는 반대로 '일단 가야 할 위치(p)를 정해두고, 제약 조건에 맞춰 위치를 수정한 뒤 속도를 역산'합니다. 이 덕분에 오차가 누적되어 물체가 멀리 날아가 버리는 현상이 근본적으로 억제되며, 매우 안정적인 시뮬레이션이 가능해집니다.

3-2. The Solver

솔버(solver)의 입력 데이터는 $M + M_{coll}$개의 제약 조건(constraints)과 각 점의 새로운 위치 예측치(estimates)인 $\mathbf{p}_1, \dots, \mathbf{p}_N$입니다. 솔버는 이러한 예측치들이 모든 제약 조건을 만족하도록 수정하는 역할을 합니다. 이때 발생하는 방정식 시스템(system of equations)은 비선형(non-linear)적인 특성을 갖는데, 아주 단순한 거리 제약 조건(distance constraint)인 $C(\mathbf{p}_1, \mathbf{p}_2) = |\mathbf{p}_1 − \mathbf{p}_2| − d$조차 비선형 방정식을 생성할 정도입니다. 또한 부등식(inequality) 타입의 제약 조건들은 결과적으로 부등식 체계를 형성합니다.

이러한 일반적인 방정식 및 부등식 집합을 해결하기 위해 본 논문에서는 가우스-자이델 방식의 반복(Gauss-Seidel-type iteration) 기법을 사용합니다. 본래 가우스-자이델 알고리즘(GS)은 선형 시스템(linear system) 전용이지만, 우리는 각 제약 조건을 하나씩 순차적으로 독립 연구하여 해결한다는 핵심 아이디어를 차용했습니다. 다만 GS와 다른 점은 개별 제약 조건을 해결하는 과정 자체가 비선형 연산(non linear operation)이라는 것입니다. 우리는 모든 제약 조건을 반복적으로 순회하며, 오직 해당 제약 조건만을 기준으로 입자들을 유효한 위치로 투영(project)합니다.

이 방식은 야코비 방식의 반복(Jacobi-type iteration)과 달리 위치 수정 결과가 다음 계산에 즉시 반영된다는 장점이 있습니다. 덕분에 압력파(pressure waves)가 단일 솔버 단계 내에서 물체 전체로 전파될 수 있어 수렴(convergence) 속도가 획기적으로 빨라집니다. 다만 이러한 효과는 제약 조건을 해결하는 순서(order)에 영향을 받으며, 제약 조건이 과도하게 설정된(over-constrained) 상황에서 이 순서가 일정하게 유지되지 않을 경우 시스템에 진동(oscillations)이 발생할 수 있습니다.

💡 가우스-자이델(Gauss-Seidel) 방식의 이점
솔버가 제약 조건을 해결할 때, '한꺼번에' 계산하지 않고 '하나씩 순서대로' 바로잡습니다. 예를 들어, 1번 점을 옮기면 그 결과가 즉시 2번 점 계산에 반영됩니다. 이는 마치 실타래를 한쪽에서 당기면 반대편까지 즉각 장력이 전달되는 것과 비슷하여, 시뮬레이션이 매우 빠르게 안정적인 상태(수렴)에 도달하게 도와줍니다.

3-3. Constraint Projection

제약 조건(constraint)에 따라 일련의 점들을 투영(projecting)한다는 것은 해당 제약 조건을 만족하도록 점들을 이동시키는 것을 의미합니다. 시뮬레이션 루프(simulation loop) 내부에서 점들을 직접 이동시키는 것과 관련하여 가장 중요한 이슈는 선운동량(linear momentum)과 각운동량(angular momentum)의 보존(conservation)입니다.

$\Delta \mathbf{p}_i$를 투영에 의한 정점(vertex) $i$의 변위(displacement)라고 할 때, 만약

• (1) $\sum_i m_i \Delta \mathbf{p}_i = 0$

이라면 선운동량이 보존되며 이는 질량 중심(center of mass)을 보존하는 것과 같습니다. 또한

• (2) $\sum_i \mathbf{r}_i \times m_i \Delta \mathbf{p}_i = 0$

라면 각운동량이 보존되며, 여기서 $\mathbf{r}_i$는 임의의 공통 회전 중심에 대한 $\mathbf{p}_i$의 거리입니다.

만약 투영이 이러한 제약 조건 중 하나를 위반하면, 객체를 끌어당기거나 회전시키는 외부 힘(external forces)처럼 작용하는 소위 유령 힘(ghost forces)이 도입됩니다. 하지만 오직 내부 제약 조건(internal constraints)들만이 운동량(momenta)을 보존할 필요가 있으며, 충돌(collision)이나 부착(attachment) 제약 조건은 객체에 전역적인 영향(global effects)을 미치는 것이 허용됩니다.

💡 '유령 힘(Ghost Forces)'을 막는 수학적 원리
물체를 시뮬레이션할 때 우리가 의도하지 않은 힘이 생겨서 물체가 스스로 회전하거나 공중으로 떠오르면 안 됩니다. 식 (5)와 (9)에서 그라디언트($\nabla C$) 방향으로만 점을 옮기는 이유는, 그 방향이 물체의 전체적인 회전이나 이동을 건드리지 않으면서 오직 '제약 조건(길이, 각도 등)'만 만족시키는 최단 거리이기 때문입니다.

우리가 제약 조건 투영(constraint projection)을 위해 제안하는 방법은 내부 제약 조건(internal constraints)에 대해 두 운동량(momenta)을 모두 보존합니다. 다시 한번 말하지만, 힘 기반 방법(force based methods)들이 에너지 항(energy term)을 통해 힘을 유도하는 것과 달리, 점 기반 접근 방식(point based approach)은 제약 함수(constraint function)를 직접 사용할 수 있다는 점에서 더 직접적입니다.
제약 함수 $C$와 강성(stiffness) $k$를 가지며 점 $\mathbf{p}_1, \dots, \mathbf{p}_n$에 작용하는 카디널리티(cardinality) $n$인 제약 조건을 살펴봅시다. 우리는 $\mathbf{p}$를 $[\mathbf{p}_1^T, \dots, \mathbf{p}_n^T]^T$로 결합된(concatenation) 벡터라고 정의합니다. 내부 제약 조건(internal constraints)의 경우, $C$는 강체 모드(rigid body modes), 즉 병진(translation) 및 회전(rotation)과 독립적입니다. 이는 점들을 회전시키거나 병진 이동시켜도 제약 함수의 값이 변하지 않음을 의미합니다. 따라서 그라디언트(gradient) $\nabla_{\mathbf{p}}C$는 최대 변화의 방향이기 때문에 강체 모드(rigid body modes)에 수직(perpendicular)입니다. 만약 보정값(correction) $\Delta\mathbf{p}$가 $\nabla_{\mathbf{p}}C$를 따라 선택된다면, 모든 질량이 동일할 시 두 운동량은 자동으로 보존됩니다(질량이 다른 경우는 나중에 다룹니다). 주어진 $\mathbf{p}$에 대해 우리는 $C(\mathbf{p}+\Delta\mathbf{p}) = 0$을 만족하는 보정값 $\Delta\mathbf{p}$를 찾고자 하며, 이 방정식은

• (3) $C(\mathbf{p}+\Delta\mathbf{p}) \approx C(\mathbf{p}) + \nabla_{\mathbf{p}}C(\mathbf{p}) \cdot \Delta\mathbf{p} = 0$

으로 근사될 수 있습니다.

$\Delta\mathbf{p}$가 $\nabla_{\mathbf{p}}C$의 방향에 있도록 제한하는 것은

• (4) $\Delta\mathbf{p} = \lambda\nabla_{\mathbf{p}}C(\mathbf{p})$

를 만족하는 스칼라(scalar) $\lambda$를 선택하는 것을 의미합니다.

식 (4)를 식 (3)에 대입하여 $\lambda$에 대해 풀고,
이를 다시 식 (4)에 대입하면 $\Delta\mathbf{p}$에 대한 최종 공식인

• (5) $\Delta\mathbf{p} = -\frac{C(\mathbf{p})}{|\nabla_{\mathbf{p}}C(\mathbf{p})|^2} \nabla_{\mathbf{p}}C(\mathbf{p}) \quad$

가 도출됩니다.
이는 단일 제약 조건(single constraint)에 의해 주어진 비선형 방정식(non-linear equation)의 반복적 해법(iterative solution)을 위한 일반적인 뉴턴-랩슨 단계(Newton-Raphson step)입니다.
개별 점(individual point) $\mathbf{p}_i$의 보정값(correction)에 대해서는

• (6) $\Delta\mathbf{p}_i = -s \nabla_{\mathbf{p}_i}C(\mathbf{p}_1, \dots, \mathbf{p}_n) \quad$

의 식을 얻게 됩니다.

여기서 스케일링 인자(scaling factor)

• (7) $s = \frac{C(\mathbf{p}_1, \dots, \mathbf{p}n)}{\sum_j |\nabla{\mathbf{p}_j}C(\mathbf{p}_1, \dots, \mathbf{p}_n)|^2} \quad$

는 모든 점에 대해 동일하게 적용됩니다.
만약 점들이 개별적인 질량을 가진다면, 우리는 보정값(corrections) $\Delta\mathbf{p}_i$를 역질량(inverse masses) $w_i = 1/m_i$에 따라 가중치를 부여합니다.
이 경우, 예상하는 바와 같이 무한한 질량(infinite mass), 즉 $w_i = 0$인 점은 움직이지 않게 됩니다.

이제 식 (4)는 $\Delta\mathbf{p}i = \lambda w_i \nabla{\mathbf{p}_i}C(\mathbf{p})$로 대체되며, 이는 스케일링 인자(scaling factor)에 대해

• (8) $s = \frac{C(\mathbf{p}_1, \dots, \mathbf{p}n)}{\sum_j w_j |\nabla{\mathbf{p}_j}C(\mathbf{p}_1, \dots, \mathbf{p}_n)|^2} \quad$

를 산출하고, 최종 보정값(final correction)으로

• (9) $\Delta\mathbf{p}i = -s w_i \nabla{\mathbf{p}_i}C(\mathbf{p}_1, \dots, \mathbf{p}_n) \quad$

를 도출하게 됩니다.

예를 들어, 거리 제약 함수(distance constraint function) $C(\mathbf{p}_1, \mathbf{p}_2) = |\mathbf{p}_1 - \mathbf{p}_2| - d$를 고려해 봅시다.
각 점에 대한 도함수(derivative)는 $\mathbf{n} = \frac{\mathbf{p}_1 - \mathbf{p}_2}{|\mathbf{p}_1 - \mathbf{p}_2|}$라고 할 때, $\nabla_{\mathbf{p}_1}C(\mathbf{p}_1, \mathbf{p}_2) = \mathbf{n}$ 및 $\nabla_{\mathbf{p}_2}C(\mathbf{p}_1, \mathbf{p}_2) = -\mathbf{n}$이 됩니다. 따라서 스케일링 인자(scaling factor) $s$는 $s = \frac{|\mathbf{p}_1 - \mathbf{p}_2| - d}{w_1 + w_2}$로 계산되며, 최종 보정값(final corrections)은

• (10) $\Delta\mathbf{p}_1 = -\frac{w_1}{w_1 + w_2}(|\mathbf{p}_1 - \mathbf{p}_2| - d)\frac{\mathbf{p}_1 - \mathbf{p}_2}{|\mathbf{p}_1 - \mathbf{p}_2|} \quad$

• (11) $\Delta\mathbf{p}_2 = +\frac{w_2}{w_1 + w_2}(|\mathbf{p}_1 - \mathbf{p}_2| - d)\frac{\mathbf{p}_1 - \mathbf{p}_2}{|\mathbf{p}_1 - \mathbf{p}_2|} \quad$

과 같이 도출됩니다.

그림 2 - 제약 조건 $C(\mathbf{p}_1, \mathbf{p}_2) = |\mathbf{p}_1 - \mathbf{p}_2| - d$의 투영(Projection)을 보여줍니다. 보정값(corrections) $\Delta \mathbf{p}_i$는 역질량(inverse masses) $w_i = 1/m_i$에 따라 가중치가 부여되어 계산됩니다.

이것들은 거리 제약 조건(distance constraints)의 투영(projection)을 위해 [Jak01]에서 제안된 공식들입니다 (그림 2 참조).
이 공식들은 본 논문에서 제시한 일반적인 제약 조건 투영 방법(general constraint projection method)의 특수한 사례(special case)로서 도출됩니다.

우리는 지금까지 제약 조건의 타입(type)과 강성(stiffness) $k$를 고려하지 않았습니다. 타입 처리(Type handling)는 간단합니다. 만약 타입이 등식(equality)이라면 우리는 항상 투영(projection)을 수행하고, 만약 타입이 부등식(inequality)이라면 $C(\mathbf{p}_1, \dots, \mathbf{p}_n) < 0$인 경우에만 투영이 수행됩니다.
강성 파라미터(stiffness parameter)를 포함시키는 데는 몇 가지 방법이 있는데, 가장 단순한 변형은 보정값(corrections) $\Delta\mathbf{p}$에 $k \in [0 \dots 1]$를 곱하는 것입니다. 하지만 솔버의 다중 반복 루프(multiple iteration loops)에 대해 $k$의 효과는 비선형(non-linear)적입니다. $n_s$번의 솔버 반복 후 단일 거리 제약 조건에 대해 남는 오차(remaining error)는 $\Delta\mathbf{p}(1 - k)^{n_s}$가 됩니다. 선형 관계(linear relationship)를 얻기 위해 우리는 보정값에 $k$를 직접 곱하는 대신 $k' = 1 - (1 - k)^{1/n_s}$를 곱합니다. 이 변환을 통해 오차는 $\Delta\mathbf{p}(1 - k')^{n_s} = \Delta\mathbf{p}(1 - k)$가 되며, 따라서 원하는 대로 $k$에 선형적으로 의존하고 $n_s$에는 독립적이게 됩니다.
하지만 결과로 나타나는 재질 강성(material stiffness)은 여전히 시뮬레이션의 타임 스텝(time step)에 의존합니다. 실시간 환경(Real time environments)은 보통 고정 타임 스텝(fixed time steps)을 사용하며, 이 경우 이러한 의존성은 문제가 되지 않습니다.

3-4. Collision Detection and Response

위치 기반 접근 방식(position based approach)의 한 가지 장점은 충돌 응답(collision response)이 얼마나 단순하게 실현될 수 있는지에 있습니다. 시뮬레이션 알고리즘의 (8)번 라인에서 $M_{coll}$개의 충돌 제약 조건(collision constraints)이 생성됩니다. 객체 표현에 의해 주어진 처음 $M$개의 제약 조건은 시뮬레이션 전반에 걸쳐 고정(fixed)되어 있는 반면, 추가적인 $M_{coll}$개의 제약 조건은 매 타임 스텝(time step)마다 처음부터 생성됩니다. 충돌 제약 조건의 수 $M_{coll}$은 가변적이며 충돌하는 정점(vertices)의 수에 의존합니다. 여기서는 연속 충돌(continuous collisions)과 정적 충돌(static collisions) 모두 처리될 수 있습니다.
연속 충돌 처리(continuous collision handling)를 위해 각 정점 $i$에 대하여 레이(ray) $\mathbf{x}_i \to \mathbf{p}_i$를 테스트합니다. 만약 이 레이가 객체 내부로 진입한다면, 진입점(entry point) $\mathbf{q}_c$와 해당 위치의 표면 법선(surface normal) $\mathbf{n}_c$를 계산합니다. 이후 제약 함수 $C(\mathbf{p}) = (\mathbf{p} - \mathbf{q}_c) \cdot \mathbf{n}_c$와 강성(stiffness) $k = 1$을 가진 부등식 제약 조건(inequality constraint)을 리스트에 추가합니다. 만약 레이 $\mathbf{x}_i \to \mathbf{p}_i$가 완전히 객체 내부에 존재한다면, 이는 연속 충돌 감지(continuous collision detection)가 특정 시점에 실패했음을 의미하며, 이 경우 정적 충돌 처리(static collision handling)로 대체(fall back)합니다. 이때는 $\mathbf{p}_i$에 가장 가까운 표면 지점 $\mathbf{q}_s$와 표면 법선 $\mathbf{n}_s$를 구하여 제약 함수 $C(\mathbf{p}) = (\mathbf{p} - \mathbf{q}_s) \cdot \mathbf{n}_s$를 가진 제약 조건을 추가합니다. 충돌 제약 조건 생성은 솔버 루프(solver loop) 외부에서 수행되므로 시뮬레이션 속도가 훨씬 빨라지지만, 고정된 제약 조건 집합을 사용함에 따라 충돌을 놓치는 시나리오가 생길 수 있습니다. 그러나 경험상 이러한 아티팩트(artifacts)는 무시할 수 있는 수준입니다.

💡 위치 기반 충돌 처리의 간결함
PBD에서 충돌 처리는 복잡한 충격량 계산이 필요 없습니다. 그냥 "어? 뚫고 들어갔네? 그럼 밖으로 꺼내자"라고 위치 제약 조건(Static Collision)을 하나 추가하는 것으로 끝납니다.

마찰(Friction)과 반발(restitution)은 알고리즘의 (16)번 단계에서 충돌하는 정점들의 속도(velocities)를 조작함으로써 처리될 수 있습니다. 충돌 제약 조건이 생성된 각 정점의 속도는 충돌 법선(collision normal)에 수직인 방향으로 감쇠(dampened)되고, 충돌 법선 방향으로 반사(reflected)됩니다. 위에서 논의된 충돌 처리는 충돌 파트너에게 충격량(impulse)이 전달되지 않기 때문에 정적 객체(static objects)와의 충돌에 대해서만 정확합니다. 두 동적 충돌 객체(dynamic colliding objects)에 대한 올바른 응답은 두 객체 모두를 우리의 시뮬레이터로 시뮬레이션함으로써 달성될 수 있는데, 즉 우리 알고리즘의 입력인 $N$개의 정점과 $M$개의 제약 조건이 단순히 둘 이상의 독립적인 객체를 나타내게 하는 것입니다.

그 후, 만약 한 객체의 점 $\mathbf{q}$가 다른 객체의 삼각형 $\mathbf{p}_1, \mathbf{p}_2, \mathbf{p}_3$를 통과하여 이동한다면, 점 $\mathbf{q}$를 삼각형의 올바른 쪽에 유지시키는 제약 함수 $C(\mathbf{q}, \mathbf{p}_1, \mathbf{p}_2, \mathbf{p}_3) = \pm(\mathbf{q} - \mathbf{p}_1) \cdot [(\mathbf{p}_2 - \mathbf{p}_1) \times (\mathbf{p}_3 - \mathbf{p}_1)]$를 가진 부등식 제약 조건(inequality constraint)을 삽입합니다. 이 제약 함수는 강체 모드(rigid body modes)와 독립적이기 때문에, 선운동량과 각운동량을 정확하게 보존(conserve)할 것입니다. 충돌 감지(Collision detection)는 네 개의 정점이 레이(rays) $\mathbf{x}_i \to \mathbf{p}_i$로 표현되기 때문에 약간 더 복잡해집니다. 따라서 움직이는 점과 움직이는 삼각형 사이의 충돌이 감지되어야 합니다 (옷감 자체 충돌(cloth self collision)에 관한 섹션 참조).

3-5. Damping

시뮬레이션 알고리즘의 (6)번 라인에서, 속도(velocities)는 새로운 위치를 예측(prediction)하는 데 사용되기 전에 감쇠(dampened)됩니다.
어떤 형태의 감쇠(damping)든 사용될 수 있으며 문헌상에 이미 많은 감쇠 방법들이 제안되어 왔으나([NMK∗05] 참조), 여기서 우리는 다음과 같은 몇 가지 흥미로운 특성을 가진 새로운 방법을 제안합니다.

(1) $\mathbf{x}_{cm} = (\sum_i \mathbf{x}_i m_i) / (\sum_i m_i)$
(2) $\mathbf{v}_{cm} = (\sum_i \mathbf{v}_i m_i) / (\sum_i m_i)$
(3) $\mathbf{L} = \sum_i \mathbf{r}_i \times (m_i \mathbf{v}_i)$
(4) $\mathbf{I} = \sum_i \tilde{\mathbf{r}}_i \tilde{\mathbf{r}}_i^T m_i$
(5) $\boldsymbol{\omega} = \mathbf{I}^{-1} \mathbf{L}$
(6) forall vertices $i$
(7) $\quad$ $\Delta \mathbf{v}_i = \mathbf{v}_{cm} + \boldsymbol{\omega} \times \mathbf{r}_i - \mathbf{v}_i$
(8) $\quad$ $\mathbf{v}_i \leftarrow \mathbf{v}_i + k_{damping} \Delta \mathbf{v}_i$
(9) endfor

여기서 $\mathbf{r}_i = \mathbf{x}_i - \mathbf{x}_{cm}$이며, $\tilde{\mathbf{r}}_i$는 $\tilde{\mathbf{r}}_i \mathbf{v} = \mathbf{r}_i \times \mathbf{v}$라는 성질을 만족하는 $3 \times 3$ 행렬(matrix)입니다. 또한 $k_{damping} \in [0 \dots 1]$는 감쇠 계수(damping coefficient)를 나타냅니다.
알고리즘의 (1)~(5)번 라인에서 시스템의 전역 선속도(global linear velocity) $\mathbf{x}_{cm}$와 각속도(angular velocity) $\boldsymbol{\omega}$를 계산한 후, (6)~(9)번 라인을 통해 전역 운동(global motion) $\mathbf{v}_{cm} + \boldsymbol{\omega} \times \mathbf{r}_i$으로부터 발생하는 개별 속도 $\mathbf{v}_i$의 편차(deviations) $\Delta \mathbf{v}_i$만을 감쇠시킵니다.
결과적으로 $k_{damping} = 1$인 극단적인 상황에서는 오직 전역 운동만이 유지되어 정점(vertices) 집합이 마치 강체(rigid body)처럼 거동하게 됩니다. 임의의 $k_{damping}$ 값을 사용할 경우, 정점들의 전역 운동(global motion) 자체에는 영향을 주지 않으면서 속도들을 전역적으로 감쇠시킬 수 있습니다.

💡 감쇠(Damping) 기법
보통 감쇠를 적용하면 물체가 허공에서 멈춰버리는 느낌을 줍니다. 하지만 여기서 제안하는 방식은 '물체 전체의 이동 속도'와 '회전 속도'는 건드리지 않고, 개별 입자들이 제멋대로 떨리는 속도만 줄입니다.

3-6. Attachments

위치 기반 접근 방식(position based approach)을 사용하면, 정점(vertices)을 정적 객체(static objects)나 운동학적 객체(kinematic objects)에 부착하는 것이 매우 간단합니다.
정점의 위치를 단순히 정적 대상 위치로 설정하거나, 운동학적 객체의 위치와 일치하도록 매 타임 스텝(time step)마다 업데이트하기만 하면 됩니다.
이 정점을 포함하는 다른 제약 조건(constraints)들이 해당 정점을 이동시키지 못하도록 확실히 하기 위해, 해당 정점의 역질량(inverse mass) $w_i$를 0으로 설정합니다.

💡 무한한 질량($w=0$)의 역할
특정 점을 벽에 고정하고 싶을 때, 역질량을 0으로 만드는 것만큼 확실한 방법은 없습니다. 수식 (9)를 보면 $w=0$일 때 보정값($\Delta p$)이 무조건 0이 됩니다. 어떤 강력한 제약 조건이 당겨도 이 점은 꿈쩍도 하지 않게 되는데, 이를 이용해 캐릭터의 손에 옷감을 쥐여주거나 벽에 깃발을 고정하는 기능을 아주 쉽게 구현할 수 있습니다.

---

4. Cloth Simulation

우리는 게임을 위한 실시간 옷감 시뮬레이터(real time cloth simulator)를 구현하기 위해 점 기반 역학 프레임워크(point based dynamics framework)를 사용했습니다. 이 섹션에서 우리는 옷감 특유의 문제들(cloth specific issues)을 다룰 것이며, 이를 통해 이전 섹션에서 소개된 일반적인 개념(general concepts)들에 대한 구체적인 예시를 제공할 것입니다.

4-1. Representation of Cloth

우리의 옷감 시뮬레이터(cloth simulator)는 임의의 삼각형 메쉬(triangle meshes)를 입력으로 받습니다. 입력 메쉬에 부여하는 유일한 제한은 그것이 매니폴드(manifold)를 나타내야 한다는 것, 즉 각 에지(edge)가 최대 두 개의 삼각형에 의해 공유되어야 한다는 것입니다. 메쉬의 각 노드(node)는 시뮬레이션되는 정점(vertex)이 됩니다. 사용자는 면적당 질량($[kg/m^2]$)으로 주어지는 밀도(density) $\rho$를 제공하며, 정점의 질량은 인접한 각 삼각형 질량의 3분의 1을 합산한 값으로 설정됩니다.각 에지에 대해, 우리는 제약 함수

• $C_{stretch}(\mathbf{p}_1, \mathbf{p}_2) = |\mathbf{p}_1 - \mathbf{p}_2| - l_0$

와 강성(stiffness) $k_{stretch}$ 및 등식(equality) 타입을 가진 인장 제약 조건(stretching constraint)을 생성합니다. 스칼라 $l_0$는 에지의 초기 길이이며 $k_{stretch}$는 사용자가 제공하는 옷감의 인장 강성(stretching stiffness)을 정의하는 전역 파라미터(global parameter)입니다.

인접한 삼각형 쌍 $(\mathbf{p}_1, \mathbf{p}_3, \mathbf{p}2)$와 $(\mathbf{p}_1, \mathbf{p}_2, \mathbf{p}_4)$에 대해, 우리는 다음과 같은 제약 함수(constraint function)

• $C_{bend}(\mathbf{p}_1, \mathbf{p}_2, \mathbf{p}_3, \mathbf{p}_4) = \text{acos} \left( \frac{(\mathbf{p}_2 - \mathbf{p}_1) \times (\mathbf{p}_3 - \mathbf{p}_1)}{|(\mathbf{p}_2 - \mathbf{p}_1) \times (\mathbf{p}_3 - \mathbf{p}_1)|} \cdot \frac{(\mathbf{p}_2 - \mathbf{p}_1) \times (\mathbf{p}_4 - \mathbf{p}_1)}{|(\mathbf{p}_2 - \mathbf{p}_1) \times (\mathbf{p}_4 - \mathbf{p}_1)|} \right) - \phi_0$

와 강성(stiffness) $k_{bend}$ 및 등식(equality) 타입을 가진 굽힘 제약 조건(bending constraint)을 생성합니다. 스칼라 $\phi_0$는 두 삼각형 사이의 초기 이면각(dihedral angle)이며, $k{bend}$는 옷감의 굽힘 강성(bending stiffness)을 정의하는 전역 사용자 파라미터(global user parameter)입니다 (그림 4 참조).

그림 4: 굽힘 저항(bending resistance)을 위해 제약 함수 $C(\mathbf{p}_1, \mathbf{p}_2, \mathbf{p}_3, \mathbf{p}_4) = \arccos(\mathbf{n}_1 \cdot \mathbf{n}_2) - \phi_0$가 사용됩니다. 실제 이면각(dihedral angle) $\phi$는 두 삼각형의 법선(normals) 사이의 각도로 측정됩니다.

점 $\mathbf{p}_3$와 $\mathbf{p}_4$ 사이에 단순 거리 제약 조건을 추가하는 방식보다 이 굽힘 항이 갖는 장점은 인장(stretching)과 독립적이라는 점입니다.
이는 해당 항이 에지 길이(edge lengths)와 무관하기 때문이며, 덕분에 사용자는 인장 강성은 낮지만 굽힘 저항(bending resistance)은 높은 옷감을 독립적으로 설정할 수 있습니다 (그림 3 참조).

그림 3: 우리가 제안하는 굽힘 항(bending term)을 사용하면, 굽힘(bending)과 인장(stretching)은 독립적인 파라미터(independent parameters)가 됩니다.
윗줄은 $(k_{stretching}, k_{bending}) = (1,1)$, $(\frac{1}{2}, 1)$, 그리고 $(\frac{1}{100}, 1)$인 경우를 보여줍니다.
아랫줄은 $(k_{stretching}, k_{bending}) = (1,0)$, $(\frac{1}{2}, 0)$, 그리고 $(\frac{1}{100}, 0)$인 경우를 보여줍니다.

식 (10)과 (11)은 인장 제약 조건(stretching constraints)에 대한 투영(projection)을 정의합니다.
부록 A(appendix A)에서 우리는 굽힘 제약 조건(bending constraints)을 투영하기 위한 공식들을 도출합니다.

4-2. Collision with Rigid Bodies

강체(rigid bodies)와의 충돌 처리를 위해 우리는 섹션 3.4에서 설명한 대로 진행합니다.
양방향 상호작용(two-way interactions)을 얻기 위해, 정점 $i$가 해당 물체와의 충돌로 인해 투영(projected)될 때마다 접촉 지점의 강체에 충격량(impulse) $m_i \Delta \mathbf{p}_i / \Delta t$를 가합니다. 옷감의 정점들에 대해서만 충돌 테스트를 수행하는 것은 작은 강체가 큰 옷감 삼각형(cloth triangles) 사이로 빠져나갈 수 있기 때문에 충분하지 않습니다. 따라서 강체의 볼록한 코너(convex corners)와 옷감 삼각형 사이의 충돌도 함께 테스트합니다.

4-3. Self Collision

삼각형들이 모두 거의 같은 크기를 갖는다고 가정하고, 우리는 정점-삼각형 충돌(vertex triangle collisions)을 찾기 위해 공간 해싱(spatial hashing)을 사용합니다 [THM∗03]. 만약 정점 $\mathbf{q}$가 삼각형 $\mathbf{p}_1, \mathbf{p}_2, \mathbf{p}_3$를 통과하여 이동하면, 우리는 다음과 같은 제약 함수를 사용합니다.

• (12) $C(\mathbf{q},\mathbf{p}_1,\mathbf{p}_2,\mathbf{p}_3) = (\mathbf{q}-\mathbf{p}_1)\cdot \frac{(\mathbf{p}_2 -\mathbf{p}_1)\times(\mathbf{p}_3 -\mathbf{p}_1)}{|(\mathbf{p}_2 -\mathbf{p}_1)\times(\mathbf{p}_3 -\mathbf{p}_1)|} -h, \quad$

여기서 $h$는 옷감의 두께(cloth thickness)를 나타냅니다 (그림 5 참조).

그림 5: 제약 함수 $C(\mathbf{q}, \mathbf{p}_1, \mathbf{p}_2, \mathbf{p}_3) = (\mathbf{q} - \mathbf{p}_1) \cdot \mathbf{n} - h$는 정점 $\mathbf{q}$가 삼각형 $\mathbf{p}_1, \mathbf{p}_2, \mathbf{p}_3$의 위쪽에 옷감 두께(cloth thickness) $h$만큼의 거리를 유지하며 머물도록 보장합니다.

만약 정점이 삼각형 법선(normal) 기준으로 아래쪽에서 진입한다면, 정점을 원래 있던 쪽에 유지시키기 위해 제약 함수는 다음과 같이 정의되어야 합니다.

• (13) $C(\mathbf{q},\mathbf{p}_1,\mathbf{p}_2,\mathbf{p}_3) = (\mathbf{q}-\mathbf{p}_1)\cdot \frac{(\mathbf{p}_3 -\mathbf{p}_1)\times(\mathbf{p}_2 -\mathbf{p}_1)}{|(\mathbf{p}_3 -\mathbf{p}_1)\times(\mathbf{p}_2 -\mathbf{p}_1)|} -h \quad$

이러한 제약 조건들을 투영(Projecting)하는 것은 선운동량과 각운동량을 보존하며, 이는 자체 충돌이 내부적인 과정(internal process)이기 때문에 옷감 자체 충돌(cloth self collision)에 있어 매우 필수적입니다. 그림 6은 자체 충돌이 포함된 옷감 조각의 정지 상태를 보여줍니다. 다만 옷감이 복잡하게 엉킨 상태(tangled state)가 될 경우 단순히 연속 충돌(continuous collisions)을 테스트하는 것만으로는 불충분할 수 있으며, 이 경우 [BWK03]에서 제안된 것과 같은 고급 기법들이 적용되어야 합니다.

그림 6: 이 접혀 있는 설정(folded configuration)은 안정적인 자체 충돌(self collision)과 그에 따른 반응(response)을 보여줍니다.

4-4. Cloth Balloons

닫힌 삼각형 메쉬(closed triangle meshes)의 경우, 메쉬 내부의 과압(overpressure)을 쉽게 모델링할 수 있습니다 (그림 7 참조).

그림 7: 캐릭터 내부의 과압(overpressure) 시뮬레이션.

우리는 제약 조건 집합에 강성(stiffness) $k = 1$을 적용하며, 메쉬의 모든 $N$개 정점(vertices)에 관계된 다음과 같은 제약 함수(constraint function)를 가진 등식 제약 조건(equality constraint)을 추가합니다.

• (14) $C(\mathbf{p}_1, \dots, \mathbf{p}_N) = \left( \sum_{i=1}^{n_{triangles}} (\mathbf{p}_{t_1^i} \times \mathbf{p}_{t_2^i}) \cdot \mathbf{p}_{t_3^i} \right) - k_{pressure}V_0 \quad$

여기서 $t_1^i, t_2^i, t_3^i$는 삼각형 $i$에 속한 정점들의 세 인덱스(indices)입니다. 이 합산식은 닫힌 메쉬의 실제 부피(actual volume)를 계산하는데, 이 값을 초기 부피 $V_0$에 과압 계수(overpressure factor) $k_{pressure}$를 곱한 목표 부피와 비교합니다. 이 제약 함수는 다음과 같은 그라디언트(gradients) 식을 도출합니다.

• (15) $\nabla_{\mathbf{p}_i}C = \sum_{j:t_1^j=i} (\mathbf{p}_{t_2^j} \times \mathbf{p}_{t_3^j}) + \sum_{j:t_2^j=i} (\mathbf{p}_{t_3^j} \times \mathbf{p}_{t_1^j}) + \sum_{j:t_3^j=i} (\mathbf{p}_{t_1^j} \times \mathbf{p}_{t_2^j}) \quad$

최종 투영 오프셋(projection offsets) $\Delta \mathbf{p}_i$를 얻으려면, 이 그라디언트들을 식 (7)의 스케일링 계수로 조정하고 식 (9)에 따라 질량(masses)에 따른 가중치를 부여해야 합니다.

---

5. Results

우리는 우리의 방법을 물리 시뮬레이션을 위한 게임과 같은 환경(gamelike environment)인 Rocket [Rat04]에 통합하였습니다. 제안된 방법의 특성과 성능(performance)을 분석하기 위해 다양한 실험이 수행되었으며, 이 섹션에서 제시되는 모든 테스트 시나리오는 PC Pentium 4, 3 GHz 사양에서 수행되었습니다.

• 독립적인 굽힘 및 인장(Independent Bending and Stretching)

우리가 제안하는 굽힘 항(bending term)은 에지 길이(edge lengths)가 아닌 인접한 삼각형들의 이면각(dihedral angle)에만 의존합니다. 따라서 굽힘 및 인장 저항을 독립적인 파라미터로 선택할 수 있습니다. 그림 3은 다양한 인장 강성(stretching stiffnesses)을 가진 옷감 가방을 보여주는데, 먼저 굽힘 저항을 활성화한 경우와 비활성화한 경우를 대조합니다. 윗줄의 결과에서 확인되듯, 굽힘 특성은 인장 저항에 영향을 미치지 않습니다.

• 양방향 상호작용 부착(Attachments with Two Way Interaction)

본 프레임워크는 단방향 및 양방향으로 결합된 부착 제약 조건(attachment constraints)을 모두 시뮬레이션할 수 있습니다. 그림 8에서 옷감 띠들은 상단의 정적 강체(static rigid bodies)에 단방향 제약 조건으로 고정되어 있으며, 동시에 띠들과 하단의 강체들 사이에는 양방향 상호작용(two way interaction)이 활성화되어 있습니다. 이러한 설정은 띠들의 사실적인 흔들림과 비틀림 동작을 구현합니다. 해당 장면은 6개의 강체와 3개의 옷감 조각을 포함하고 있으며, 380fps 이상의 매우 높은 성능으로 시뮬레이션 및 렌더링됩니다.

그림 8: 옷감 띠(cloth stripes)들은 상단에서 단방향 상호작용(one way interaction)을 통해 정적 강체(static rigid bodies)에 부착되어 있으며, 하단에서는 양방향 제약 조건(two way constraints)을 통해 강체들에 부착되어 있습니다.

• 실시간 자체 충돌(Real Time Self Collision)

그림 6의 옷감 조각은 1364개의 정점(vertices)과 2562개의 삼각형(triangles)으로 구성됩니다. 이 시뮬레이션은 자체 충돌 감지(self collision detection), 충돌 처리 및 렌더링을 모두 포함하여 평균 30fps로 실행됩니다. 마찰(friction) 효과는 그림 9에서 회전하는 통 내부에서 굴러가는 옷감을 통해 확인할 수 있습니다.

• 찢어짐과 안정성(Tearing and stability)

그림 10에서는 4264개의 정점과 8262개의 삼각형을 가진 옷감이 큐브에 의해 벌어지고 공에 의해 찢어지는 장면을 보여주며, 평균 47fps의 성능을 기록했습니다. 찢어짐(Tearing)은 에지의 늘어남(stretching)이 임계값(threshold)을 넘을 때 정점을 지나 에지에 수직인 분할 평면(split plane)을 생성하여 정점을 분할(split)하는 간단한 방식으로 시뮬레이션됩니다. 우리의 방법은 그림 1처럼 부풀려진 캐릭터가 회전 기어 사이로 끼어 들어가는 극단적인 상황에서도 다수의 제약 조건(constraints)과 충돌을 안정적으로 처리해냅니다.

그림 10: 한 조각의 옷감(piece of cloth)이 부착된 큐브(attached cube)에 의해 찢어지고(torn open), 던져진 공(thrown ball)에 의해 갈기갈기 찢겨 나갑니다(ripped apart).

• 복잡한 시뮬레이션 시나리오(Complex Simulation Scenarios)

이 방법은 복잡한 시뮬레이션 환경에 특히 적합합니다(그림 12 참조). 애니메이션 캐릭터나 복잡한 지형의 게임 레벨과의 상호작용이 많더라도, 여러 개의 옷감 시뮬레이션을 실시간 대화형 속도(interactive speed)로 처리할 수 있습니다.

그림 12: 옷감 조각(pieces of cloth)들과 애니메이션되는 게임 캐릭터(animated game character, 왼쪽), 기하학적으로 복잡한 게임 레벨(geometrically complex game level, 가운데), 그리고 수백 개의 시뮬레이션된 식물 잎사귀(plant leaves, 오른쪽) 사이의 광범위한 상호작용(extensive interaction).

---

6. Conclusions

우리는 제약 함수(constraint functions)를 통해 공식화된 일반적인 제약 조건(general constraints)들을 처리할 수 있는 위치 기반 역학 프레임워크(position based dynamics framework)를 제시하였습니다. 위치 기반 접근 방식(position based approach)을 사용하면 시뮬레이션 중에 객체를 직접 조작하는 것이 가능합니다. 이는 충돌(collisions), 부착 제약 조건(attachment constraints) 및 명시적 적분(explicit integration)의 처리를 크게 단순화하며, 애니메이션 장면(animated scene)을 직접적이고 즉각적으로 제어할 수 있게 해줍니다.

우리는 이 프레임워크를 기반으로 강체와의 양방향 상호작용(two way interaction), 옷감의 자체 충돌 및 반응(self collision and response), 그리고 동적 강체(dynamic rigid bodies)에 옷감 조각을 부착하는 기능 등을 제공하는 강건한 옷감 시뮬레이터(robust cloth simulator)를 구현하였습니다.

---

7. Future Work

이 논문에서 다루지 않은 주제는 강체 시뮬레이션(rigid body simulation)입니다. 하지만 우리가 제시한 접근 방식은 강체 객체들도 처리할 수 있도록 매우 쉽게 확장될 수 있습니다. 일반적인 강체 솔버들이 충돌 해결을 위해 선운동량 및 각운동량 충격량(impulses) 세트를 계산하는 것과 달리, (PBD 방식에서는) 접촉 지점에서 물체에 직접적인 이동(movements)과 회전(rotations)을 적용한 뒤, 솔버가 완료된 후 그에 맞춰 선속도와 각속도를 조정하는 방식을 취하게 될 것입니다.