Modelarea și rezolvarea constrângerilor fizice

A constrângere este o restricție asupra unui obiect modelat fizic în cadrul unei simulări. În general, un obiect începe cu șase grade de libertate, reprezentând capacitatea sa de a se deplasa și de a se roti în lumea simulată; prin limitarea acestor grade de libertate în diferite moduri, putem obține multe efecte interesante și atrăgătoare.

Pe măsură ce CPU-ul computerelor moderne devine din ce în ce mai puternic, se pot face calcule pentru modelarea și rezolvarea multor scenarii fizice interesante. Constrângerile sunt o modalitate generalizată și matematică de a produce astfel de scenarii. Cu toții putem mulțumi lui Erin Catto pentru lucrarea sa inițială pe această temă: Dinamica iterativă cu coerență temporală. Voi folosi această lucrare ca referință la scrierea acestui post, așa că vă sugerez cel puțin să faceți o privire înainte de a citi, chiar dacă numai din respect pentru munca lui Erin (și a lucrărilor sale de resurse enumerate și mulțumiri speciale de către fiecare dintre contribuabili lor ).


Glosar

Există câțiva termeni diferiți folosiți în mod frecvent atunci când vorbim despre motoarele fizicii. Pentru claritate, iată un scurt glosar pe care îl puteți folosi ca referință la citirea acestui post:

  • Corp rigid: Un obiect simulat în cadrul unui motor de fizică. Obiectul se presupune a fi complet rigid, ca un diamant. De cele mai multe ori, corpurile rigide se învârt în jurul lor, alunecă unul pe celălalt și alcătuiesc lumea fizică a unui joc.
  • Grad de libertate: Una dintre cele șase valori scalare utilizate pentru a reprezenta poziția și orientarea corpului rigid în lume. Trei valori scalare reprezintă poziția și trei reprezintă orientarea.
  • Constrângere: O limitare pusă pe unul sau mai multe grade de libertate a unui corp rigid. Aproape toate constrângerile sunt pereche, ceea ce înseamnă că ele afectează două corpuri rigide în tandem.
  • Comun: O articulație este o constrângere pereche, deși tipul de constrângere este nu o constrângere de interpenetare. Toate celelalte tipuri de constrângeri sunt denumite în comun.
  • Contact Constrângere: O constrângere pereche care nu este o articulație, adică o constrângere de interpenetare.
  • Contactați-l pe Normal: Direcția în care trebuie rezolvată o constrângere de interpenetare. Aceasta este de obicei determinată de detectarea coliziunilor.

Cerințe preliminare

Sunt necesare câteva premise pentru a profita din plin de acest articol. Cu toate acestea, cititorul general se poate bucura încă mult de materiale, deși este mai bine să aibă o înțelegere de bază a:

  • Algebră
  • Calcul
  • Calculul vectorilor
  • Intermediar la abilități avansate de programare (pentru a scrie codul)

Tipuri de constrângeri

Deoarece o constrângere limitează gradul de libertate, să vedem ce face corpul rigid atunci când folosește șase dintr-o dată:


Corpul rigid de mai sus utilizează trei grade de libertate de a se traduce prin lume. Ultimele trei sunt utilizate pentru a schimba în mod constant orientarea tuturor celor trei axe de rotație.

Să aruncăm o privire la câteva exemple diferite de constrângeri. Cea mai cunoscută constrângere ar fi aceea care împiedică penetrarea a două corpuri rigide. Acest tip de constrângere este activă numai atunci când două corpuri se penetrează unul pe altul și conduc cele două corpuri în afară. Odată ce această constrângere de interpenetare este activă, este ușor de observat că gradele de libertate ale corpurilor rigide devin limitate și afectate în așa fel încât să producă un rezultat interesant (rezultatul interesant fiind faptul că cele două obiecte se pot ciocni unul cu celălalt ):


Lumea Bună a constrângerilor ar fi distanta de la distanta, unde două puncte pe două corpuri rigide sunt constrânse să fie o distanță exactă unul față de celălalt. Vă puteți imagina o tijă fără masă care unește două puncte împreună, unde această tijă nu poate fi întinsă sau comprimată:


Există multe tipuri de constrângeri pentru tot felul de comportamente interesante, printre care: frecare, prismatică, revoltă, sudură, unghi și multe altele.


Înțelegerea ecuațiilor de constrângere

În formă generală, o constrângere este o ecuație scalară egală cu o anumită valoare (de obicei zero).

\ Begin ecuație
C (l_1, a_1, l2, a_2) = 0
\ Label EQ1
\ End ecuație

Termenii \ (l \) și \ (a \) din \ eqref eq1 sunt propria mea notație: \ (l \) se referă la liniar în timp ce \ (a \) se referă la unghiular. Indiciile 1 și 2 se referă la cele două obiecte din cadrul constrângerii. După cum puteți vedea, există intrări liniare și unghiulare la o ecuație de constrângere și fiecare trebuie să fie o valoare scalară.

Să facem un pas înapoi să ne uităm la constrângerea de distanță. Constrângerea de distanță vrea să conducă distanța dintre două puncte de ancorare pe două corpuri pentru a fi egală cu o valoare scalară:

\ Begin ecuație
C (l_1, a_1, l2, a_2) = \ frac 1 2 [(P_2 - P_1) ^ 2 - L ^ 2]
\ Label eq2
\ End ecuație

\ (L \) este lungimea tijei care leagă ambele corpuri; \ (P_1 \) și \ (P_2 \) sunt pozițiile celor două corpuri.

În forma sa actuală, această constrângere este o ecuație de poziție. Acest tip de ecuație de poziționare este neliniar, ceea ce face ca rezolvarea să fie foarte dificilă. O metodă de rezolvare a acestei ecuații poate fi în loc de a determina constrângerea de poziție (în funcție de timp) și de a folosi o constrângere a vitezei. Rezultatele ecuațiilor de viteză sunt lineare, ceea ce le face solvabile. Soluțiile pot fi apoi integrate folosind un fel de integrator înapoi în forma pozițională.

În formă generală, o constrângere a vitezei are forma:

\ Begin ecuație
\ dot C (l_1, a_1, l_2, a_2) = 0
\ Label EQ3
\ End ecuație

Punctul deasupra lui \ (C \ in \ eqref eq3 se refera la derivatul lui \ (C \) in raport cu timpul. Aceasta este notație comună atunci când se ocupă de studiul fizicii.

În timpul derivatului, un nou termen \ (J \) apare prin regula:

\ Begin ecuație
\ dot C (l_1, a_1, l_2, a_2) = JV = 0
\ Label EQ4
\ End ecuație

Derivatul de timp al lui \ (C \) creează un vector de viteză și un Jacobian. Jacobianul este o matrice 1x6 care conține valori scalare corespunzătoare fiecărui grad de libertate. Într-o constrângere pereche, un Jacobian va conține de obicei 12 elemente (suficient pentru a conține termenii \ (l \) și \ (a \) pentru ambele corpuri \ (A \) și \ (B \).

Un sistem de constrângeri poate forma a comun. O comună poate conține numeroase constrângeri care restricționează gradele de libertate în diferite moduri. În acest caz, Jacobianul va fi o matrice în care numărul de rânduri este egal cu numărul de constrângeri active în sistem.

Jacobianul este derivat offline, de mână. Odată ce un Jacobian este dobândit, codul de calcul și de utilizare a Jacobian poate fi creat. După cum puteți vedea din \ eqref eq4, viteza \ (V \) este transformată din spațiul cartezian în spațiul de constrângere. Acest lucru este important deoarece, în spațiul de constrângere, originea este cunoscută. De fapt, orice țintă poate fi cunoscută. Acest lucru înseamnă că orice constrângere poate fi derivată pentru a da un Jacobian care poate transforma forțele din spațiul cartesian în spațiul constrângerii.

În spațiul de constrângere, având în vedere un scalar țintă, ecuația se poate deplasa spre sau departe de ținta. Soluțiile pot fi obținute cu ușurință în spațiu de constrângere pentru a muta starea curentă a unui corp rigid către o stare țintă. Aceste soluții pot fi apoi transformate din spațiul de constrângere înapoi în spațiul cartezian, cum ar fi:

\ Begin ecuație
F = \ lambda J ^ T
\ Label EQ5
\ End ecuație

\ (F \) este o forță în spațiul cartesian, unde \ (J ^ T) este inversul (transpus) iacobian. \ (\ lambda \) (lambda) este un multiplicator scalar.

Gândiți-vă la Iacobian ca vector de viteză, unde fiecare rând este un vector în sine (două valori scalare în 2D și trei valori scalare în 3D):

\ Begin ecuație
J = \ încep bmatrix
l_1 \\
a_1 \\
l_2 \\
a_2 \\
\ End bmatrix
\ Label EQ6
\ End ecuație

Multiplicarea \ (V \) cu \ (J \) ar implica matematic multiplicarea matricei. Totuși, majoritatea elementelor sunt zero și de aceea tratăm Iacobianul ca vector. Acest lucru ne permite să definim propria noastră operație de calcul \ (JV \), ca în \ eqref eq4.

\ Begin ecuație
JV = \ încep bmatrix
l_1 & a_1 & l_2 & a_2
\ End bmatrix
\ Begin bmatrix
v_1 \\
ω_1 \\
v_2 \\
ω_2 \\
\ End bmatrix
\ Label EQ7
\ End ecuație

Aici, \ (v \) reprezintă viteza liniară, iar \ (ω \) (omega) reprezintă viteza unghiulară. \ eqref eq7 pot fi scrise ca câteva produse punctuale și multiplicări pentru a oferi un calcul mai eficient comparativ cu multiplicarea completă a matricelor:

\ Begin ecuație
JV = l_1 \ cdot v_1 + a_1 \ cdot ω_1 + l_2 \ cdot v_2 + a_2 \ cdot ω_2
\ Label EQ8
\ End ecuație

Jacobianul poate fi considerat ca un vector de direcție în spațiul de constrângere. Această direcție întotdeauna indică spre țintă în direcția care necesită cea mai mică lucrare de făcut. Din moment ce această "direcție" Jacobiană este derivată offline, tot ceea ce trebuie rezolvat este magnitudinea forței care trebuie aplicată pentru a menține constrângerea. Această magnitudine se numește \ (\ lambda \). \ (\ lambda \) poate fi cunoscută sub denumirea de Lagrange Multiplier. Eu însumi nu am studiat în mod formal Mecanica Lagrangiană, totuși un studiu al Mecanicii Lagrangiane nu este necesar pentru a implementa pur și simplu constrângeri. (Sunt dovada acestui fapt!) \ (\ Lambda \) pot fi rezolvate folosind a constrângere solver (mai multe despre aceasta mai târziu).


Rezolvarea pentru Iacoveni

În lucrarea lui Erin Catto, există o schiță simplă pentru iacobienii care au derivat mâna. Pașii sunt:

  1. Începeți cu ecuația de constrângere \ (C \)
  2. Calculați timpul derivat \ (\ dot C \)
  3. Izolați toți termenii de viteză
  4. Identificați \ (J \) prin inspecție

Singura parte dificilă este calculul derivatului, iar acest lucru poate veni cu practica. În general, constrângerile derivate din mână sunt dificile, dar devin mai ușoare în timp.

Să derivăm un Jacobian valid pentru a fi folosit în rezolvarea unei constrângeri la distanță. Putem începe la pasul 1 cu \ eqref eq2. Iată câteva detalii pentru Pasul 2:

\ Begin ecuație
\ dot C = (P_2 - P_1) (\ punct P _2 - \ punct P _1)
\ Label EQ9
\ End ecuație

\ Begin ecuație
\ dot C = (P_2 - P_1) ((v_2 + ω_2 \ ori r_2) - (v_1 + ω_1 \ ori r_1))
\ Label EQ10
\ End ecuație

\ (r_1 \) și \ (r_2 \) sunt vectori de la centrul de masă la punctul de ancorare, pentru corpurile 1 și respectiv 2.

Următorul pas este de a izola termenii de viteză. Pentru a face acest lucru, vom folosi identitatea triplă a produsului scalar:

\ Begin ecuație
(P_2 - P_1) = d
\ Label EQ11
\ End ecuație

\ Begin ecuație
\ dot C = (d \ cdot v_2 + d \ cdot ω_2 \ ori r_2) - (d \ cdot v_1 + d \ cdot ω_1 \ ori r_1)
\ Label eq12
\ End ecuație

\ Begin ecuație
\ dot C = (d \ cdot v_2 + ω_2 \ cdot r_2 \ ori d) - (d \ cdot v_1 + ω_1 \ cdot r_1 \
\ Label eq13
\ End ecuație

Ultimul pas este de a identifica Jacobianul prin inspecție. Pentru a face acest lucru, toți coeficienții tuturor termenilor de viteză (\ (V \) și \ (ω \)) vor fi utilizați ca elemente Jacobian. Prin urmare:

\ Begin ecuație
J = \ begin bmatrix -d & -r_1 \ ori d & d & r_2 \ ori d \ end bmatrix
\ Label eq14
\ End ecuație

Mai mulți Jacobi

Constrângere de contact (constrângere de interpenetare), unde \ (n \) este contactul normal:

\ Begin ecuație
J = \ begin bmatrix -n & -r_1 \ timp n & n & r_2 \ timp n \ end bmatrix
\ Label eq15
\ End ecuație

Constrângerea de fricțiune (activă în timpul penetrării), unde \ (t) este o axă de frecare (2D are o axă, 3D are două):

\ Begin ecuație
J = \ begin bmatrix -t & -r_1 \ ori t & t & r_2 \ ori t \ end bmatrix
\ Label eq16
\ End ecuație


Rezolvarea constrângerilor

Acum, că avem o înțelegere a ceea ce este o constrângere, putem vorbi despre cum să le rezolvăm. Așa cum am spus mai devreme, odată ce un Iacobian este derivat de la mâini, trebuie să rezolvăm doar pentru \ (\ lambda \). Rezolvarea unei singure constrângeri în izolare este ușoară, dar rezolvarea multor constrângeri simultan este dificilă și foarte ineficientă (computațional). Aceasta reprezintă o problemă, deoarece jocurile și simulările vor dori probabil să aibă multe constrângeri active simultan.

O metodă alternativă de rezolvare simultană a tuturor constrângerilor (rezolvarea globală) ar fi rezolvarea constrângerilor iterativ. Prin rezolvarea aproximărilor soluției și prin hrănirea soluțiilor anterioare la ecuații, putem converge la soluție.

Un astfel de solver iterativ este cunoscut sub numele de impulsuri succesive, așa cum a fost numit de Erin Catto. Impulsurile secvențiale sunt foarte asemănătoare cu Gauss Seidel proiectate. Ideea este de a rezolva toate constrângerile, unul câte unul, de mai multe ori. Soluțiile se vor invalida reciproc, însă pe parcursul mai multor iterații fiecare constrângere individuală se va converge și se poate obține o soluție globală. Asta e bine! Solverii iterativi sunt rapizi.

Odată ce se ajunge la o soluție, un impuls poate fi aplicat ambelor corpuri în constrângere pentru a impune constrângerea.

Pentru a rezolva o singură constrângere, putem folosi următoarea ecuație:

\ Begin ecuație
\ lambda = \ frac - (JV + b) JM ^ - 1 J ^ T
\ Label eq17
\ End ecuație

\ (M ^ - 1 \) este masa constrângerii; \ (b \) este părtinirea (mai multe despre aceasta mai târziu).

Aceasta este o matrice care conține masa inversă și inerția inversă a ambelor corpuri rigide în constrângere. Următoarea este masa constrângerii; notează că \ (m ^ - 1 \) este masa inversă a unui corp, în timp ce \ (I ^ - 1 \) este inerția inversă a unui corp:

\ begin ecuația M ^ - 1 =
\ Begin bmatrix
m_1 ^ - 1 & 0 & 0 & 0 \\
0 & I_1 ^ - 1 & 0 & 0 \\
0 & 0 & m_2 ^ - 1 & 0 \\
0 & 0 & 0 & I_2 ​​^ - 1
\ End bmatrix
\ Label eq18
\ End ecuație

Deși \ (M ^ - 1 \) este teoretic o matrice, vă rugăm să nu modelați-o ca atare (cea mai mare parte este zero!). În schimb, fii inteligent în privința calculelor pe care le faci.

\ (JM ^ - 1 J ^ T \) este cunoscut sub numele de masă de constrângere. Acest termen este calculat o singură dată și folosit pentru a rezolva pentru \ (\ lambda \). Calculăm pentru un sistem cum ar fi:

\ Begin ecuație
JM ^ - 1 J ^ T = (l_1 \ cdot l_1) * m_1 ^ - 1 + (l_2 \ cdot l_2) * m_2 ^ a_2 * (I_2 ^ - 1 a_2)
\ Label eq19
\ End ecuație

Rețineți că trebuie să inversați \ eqref eq19 pentru a calcula \ eqref eq17.

Informațiile de mai sus sunt tot ceea ce este necesar pentru a rezolva o constrângere! O forță în spațiul cartezian poate fi rezolvată și folosită pentru a actualiza viteza unui obiect, pentru a impune o constrângere. Rețineți că \ eqref eq5:

\ Begin ecuație
F = \ lambda J ^ T \\
V_ final = V_ inițială + m ^ - 1 * F \\
∴ \\
\ Begin bmatrix
v_1 \\
ω_1 \\
v_2 \\
ω_2 \\
\ end bmatrix + = \ încep bmatrix
m_1 ^ - 1 & 0 & 0 & 0 \\
0 & I_1 ^ - 1 & 0 & 0 \\
0 & 0 & m_2 ^ - 1 & 0 \\
0 & 0 & 0 & I_2 ​​^ - 1
\ End bmatrix \ begin bmatrix
\ lambda * l_1 \\
\ lambda * a_1 \\
\ lambda * l_2 \\
\ lambda * a_2 \\
\ End bmatrix
\ Label eq20
\ End ecuație


Constrângere

Datorită linearizării ecuațiilor de poziționare neliniară, unele informații sunt pierdute. Acest lucru are ca rezultat soluții care nu satisfac destul ecuația poziției inițiale, dar satisfac ecuațiile de viteză. Această eroare este cunoscută sub numele de constrângere de constrângere. Se poate gândi la această eroare ca rezultat al aproximării liniei tangente.

Există câteva moduri diferite de a rezolva astfel de erori, toate care aproximează eroarea și aplică o anumită formă de corecție. Cel mai simplu este cunoscut sub numele de Baumgarte.

Baumgarte este un mic adaos de energie în spațiul de constrângere și explică termenul \ (b) în ecuațiile anterioare. Pentru a răspunde la părtinire, iată o versiune modificată a lui \ eqref eq4:

\ Begin ecuație
\ dot C = JV + b = 0
\ Label eq21
\ End ecuație

Pentru a calcula un termen Baumgarte și a aplica-o ca prejudecată, trebuie să inspectăm ecuația inițială de constrângere și să identificăm o metodă potrivită pentru calcularea erorii. Baumgarte este în forma:

\ Begin ecuație
JV = - \ beta C
\ Label eq22
\ End ecuație

\ (\ beta \) (termenul Baumgarte) este un factor care poate fi reglabil, fără unitate, dependent de simulare. Este, de obicei, între 0.1 și 0.3.

Pentru a calcula termenul de părtinire, să aruncăm o privire la ecuația constrângerii non-penetrare \ eqref eq15 înainte de a deriva cu privire la timp, unde \ (n \) este contactul normal:

\ Begin ecuație
C = \ begin bmatrix -x_1 & -r_1 & x_2 & r_2 \ end bmatrix \ cdot \ vec n
\ Label eq23
\ End ecuație

Ecuația de mai sus spune că eroarea scalară a lui \ (C \) este adâncimea de penetrare între două corpuri rigide.


Concluzie

Mulțumită lui Erin Catto și hârtiei sale despre rezolvarea constrângerilor, avem o cale de a rezolva constrângerile de poziție în ceea ce privește viteza. Multe comportamente interesante pot apărea din constrângeri și, sperăm, articolul va fi de folos pentru mulți oameni. Ca întotdeauna, nu ezitați să puneți întrebări sau să oferiți comentarii mai jos.

Vă rugăm să consultați Box2D Lite pentru o resursă de rezolvare a diferitelor tipuri de constrângeri 2D, precum și o privire asupra mai multor aspecte de implementare care nu sunt acoperite aici.