Modellazione e risoluzione dei vincoli fisici

UN costrizione è una restrizione su un oggetto modellato fisicamente all'interno di una simulazione. In generale, un oggetto inizia con sei gradi di libertà, rappresentando la sua capacità di muoversi e ruotare all'interno del mondo simulato; limitando questi gradi di libertà in vari modi, possiamo ottenere molti effetti interessanti e interessanti.

Mentre la CPU dei computer moderni diventa sempre più potente, i calcoli possono essere indirizzati alla modellazione e alla risoluzione di molti scenari fisici interessanti. I vincoli sono un modo generalizzato e matematicamente supportato per produrre tali scenari. Possiamo tutti ringraziare Erin Catto per il suo articolo iniziale su questo argomento: Iterative Dynamics with Temporal Coherence. Userò questo documento come riferimento quando scrivo questo post, quindi ti suggerisco di dare almeno una sbirciatina prima di leggere, anche se solo per rispetto per il lavoro di Erin (e la sua risorsa elencata funziona e ringraziamenti speciali da parte di ciascuno dei loro contributori ).


Glossario

Ci sono alcuni termini diversi comunemente usati quando si parla di motori fisici. Per chiarezza, ecco un breve glossario da utilizzare come riferimento durante la lettura di questo post:

  • Corpo rigido: Un oggetto simulato all'interno di un motore fisico. Si presume che l'oggetto sia completamente rigido, come un diamante. Il più delle volte, i corpi rigidi rimbalzano, scivolano l'uno sull'altro e costituiscono il mondo fisico di un gioco.
  • Grado di libertà: Uno dei sei valori scalari utilizzati per rappresentare la posizione e l'orientamento di un corpo rigido all'interno del mondo. Tre valori scalari rappresentano la posizione e tre rappresentano l'orientamento.
  • costrizione: Una limitazione posta su uno o più gradi di libertà di un corpo rigido. Quasi tutti i vincoli sono a coppie, nel senso che influenzano due corpi rigidi in tandem.
  • comune: Un giunto è un vincolo a coppie, sebbene sia il tipo di vincolo non un vincolo di compenetrazione. Tutti gli altri tipi di vincoli sono indicati come un giunto.
  • Contatto vincolo: Un vincolo a coppie che non è un giunto, il che significa che è un vincolo di compenetrazione.
  • Contatto normale: La direzione in cui risolvere un vincolo di compenetrazione. Questo di solito è determinato dal rilevamento delle collisioni.

Prerequisiti

Sono necessari alcuni prerequisiti per trarre il massimo vantaggio da questo articolo. Tuttavia, il lettore generale può ancora godere molto del materiale, sebbene sia meglio avere una conoscenza di base di:

  • Algebra
  • Calcolo
  • Calcolo vettoriale
  • Capacità di programmazione da intermedio ad avanzato (per scrivere il codice)

Tipi di vincoli

Dal momento che un vincolo limita i gradi di libertà, vediamo cosa fa un corpo rigido quando utilizza sei contemporaneamente:


Il corpo rigido di cui sopra usa tre gradi di libertà per tradurre se stesso attraverso il mondo. Gli ultimi tre sono utilizzati per cambiare costantemente l'orientamento di tutti e tre gli assi di rotazione.

Vediamo ora un paio di esempi diversi di quali siano effettivamente i vincoli. Il vincolo più familiare sarebbe quello che impedisce a due corpi rigidi di penetrare. Questo tipo di vincolo è attivo solo quando due corpi si penetrano l'un l'altro e allontana i due corpi. Una volta attivato questo vincolo di interpenetrazione, è facile vedere che i gradi di libertà dei corpi rigidi diventano limitati e influenzati in modo tale da produrre un risultato interessante (il risultato interessante è che i due oggetti possono scontrarsi tra loro ):


L'Hello World dei vincoli sarebbe il vincolo di distanza, dove due punti su due corpi rigidi sono vincolati per essere una distanza esatta l'uno dall'altro. Potete immaginare una verga senza massa che collega insieme due punti, dove questa asta non può essere allungata o compressa:


Esistono molti tipi di vincoli per tutti i tipi di comportamenti interessanti, tra cui: attrito, prisma, rivoluzione, saldatura, angolo e altro.


Comprendere le equazioni dei vincoli

In forma generale, un vincolo è un'equazione scalare uguale a qualche valore (solitamente zero).

\ Begin equation
C (l_1, a_1, l_2, a_2) = 0
\ Label EQ1
\ End equation

I termini \ (l \) e \ (a \) in \ eqref eq1 sono la mia notazione: \ (l \) si riferisce a lineare mentre \ (a \) si riferisce a angolare. Gli indici 1 e 2 si riferiscono ai due oggetti all'interno del vincolo. Come puoi vedere, esistono degli input lineari e angolari per un'equazione di vincoli e ognuno deve essere un valore scalare.

Facciamo un passo indietro per osservare il vincolo della distanza. Il vincolo di distanza vuole guidare la distanza tra due punti di ancoraggio su due corpi per essere uguale a qualche valore scalare:

\ Begin equation
C (l_1, a_1, l_2, a_2) = \ frac 1 2 [(P_2 - P_1) ^ 2 - L ^ 2] = 0
\ Label eq2
\ End equation

\ (L \) è la lunghezza dell'asta che collega entrambi i corpi; \ (P_1 \) e \ (P_2 \) sono le posizioni dei due corpi.

Nella sua forma attuale, questo vincolo è un'equazione di posizione. Questo tipo di equazione di posizione non è lineare, il che rende molto difficile risolverlo. Un metodo per risolvere questa equazione può invece derivare il vincolo di posizione (rispetto al tempo) e usare un vincolo di velocità. Le equazioni di velocità risultanti sono lineari, rendendole risolvibili. Le soluzioni possono quindi essere integrate utilizzando una sorta di integratore di nuovo in forma posizionale.

In forma generale, un vincolo di velocità ha la forma:

\ Begin equation
\ dot C (l_1, a_1, l_2, a_2) = 0
\ Label eq3
\ End equation

Il punto sopra \ (C \) in \ eqref eq3 si riferisce alla derivata di \ (C \) rispetto al tempo. Questa è una notazione comune quando si tratta dello studio della fisica.

Durante la derivata, un nuovo termine \ (J \) appare tramite la regola della catena:

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

La derivata temporale di \ (C \) crea un vettore di velocità e un Jacobiano. Jacobian è una matrice 1x6 contenente valori scalari corrispondenti a ciascun grado di libertà. In un vincolo pairwise, un Jacobian conterrà tipicamente 12 elementi (abbastanza per contenere i termini \ (l \) e \ (a \) per entrambi i corpi \ (A \) e \ (B \).

Un sistema di vincoli può formare un comune. Una giuntura può contenere molti vincoli che limitano i gradi di libertà in vari modi. In questo caso, Jacobian sarà una matrice in cui il numero di righe è uguale al numero di vincoli attivi nel sistema.

Il Jacobiano è derivato offline, a mano. Una volta acquisito un Jacobian, è possibile creare il codice per calcolare e utilizzare Jacobian. Come puoi vedere da \ eqref eq4, la velocità \ (V \) viene trasformata dallo spazio cartesiano allo spazio dei vincoli. Questo è importante perché nello spazio dei vincoli è nota l'origine. In realtà, qualsiasi obiettivo può essere conosciuto. Ciò significa che qualsiasi vincolo può essere derivato per produrre un Jacobiano in grado di trasformare le forze dallo spazio cartesiano allo spazio dei vincoli.

Nello spazio dei vincoli, dato uno scalare bersaglio, l'equazione può spostarsi verso o lontano dal bersaglio. Le soluzioni possono essere facilmente ottenute nello spazio dei vincoli per spostare lo stato corrente di un corpo rigido verso uno stato target. Queste soluzioni possono quindi essere trasformate dallo spazio dei vincoli nello spazio cartesiano in questo modo:

\ Begin equation
F = \ lambda J ^ T
\ Label EQ5
\ End equation

\ (F \) è una forza nello spazio cartesiano, dove \ (J ^ T \) è l'inverso (trasposto) Jacobiano. \ (\ lambda \) (lambda) è un moltiplicatore scalare.

Pensa allo Jacobiano come a un vettore di velocità, in cui ogni riga è un vettore stesso (di due valori scalari in 2D e tre valori scalari in 3D):

\ Begin equation
J = \ begin bmatrix
l_1 \\
a_1 \\
l_2 \\
a_2 \\
\ End bmatrix
\ Label eq6
\ End equation

Moltiplicare matematicamente \ (V \) per \ (J \) implicherebbe la moltiplicazione della matrice. Tuttavia, la maggior parte degli elementi è zero, ed è per questo che trattiamo il Jacobiano come un vettore. Questo ci permette di definire la nostra operazione per calcolare \ (JV \), come in \ eqref eq4.

\ Begin equation
JV = \ begin bmatrix
l_1 & a_1 & l_2 & a_2
\ End bmatrix
\ Begin bmatrix
v_1 \\
ω_1 \\
v_2 \\
ω_2 \\
\ End bmatrix
\ Label EQ7
\ End equation

Qui, \ (v \) rappresenta la velocità lineare e \ (ω \) (omega) rappresenta la velocità angolare. \ eqref eq7 può essere scritto come pochi prodotti e moltiplicazioni di punti per fornire un calcolo più efficiente rispetto alla moltiplicazione a matrice intera:

\ Begin equation
JV = l_1 \ cdot v_1 + a_1 \ cdot ω_1 + l_2 \ cdot v_2 + a_2 \ cdot ω_2
\ Label eq8
\ End equation

Il Jacobiano può essere pensato come un vettore di direzione nello spazio dei vincoli. Questa direzione punta sempre verso l'obiettivo nella direzione che richiede il minimo lavoro da fare. Poiché questa "direzione" Jacobiana è derivata offline, tutto ciò che deve essere risolto è l'entità della forza da applicare al fine di mantenere il vincolo. Questa grandezza è chiamata \ (\ lambda \). \ (\ lambda \) può essere conosciuto come il moltiplicatore di Lagrange. Io stesso non ho studiato formalmente la Meccanica Lagrangiana, tuttavia uno studio della Meccanica Lagrangiana non è necessario per implementare semplicemente dei vincoli. (Ne sono la prova!) \ (\ Lambda \) può essere risolto usando a limitatore di vincoli (più su questo più tardi).


Risolvere per Jacobians

Nell'articolo di Erin Catto esiste un semplice schema per gli Jacobiani che nascono a mano. I passaggi sono:

  1. Inizia con l'equazione dei vincoli \ (C \)
  2. Derivata del tempo di calcolo \ (\ dot C \)
  3. Isolare tutti i termini di velocità
  4. Identificare \ (J \) mediante ispezione

L'unica parte difficile è calcolare la derivata, e questo può venire con la pratica. In generale, i vincoli derivanti dalla mano sono difficili, ma diventano più facili con il tempo.

Ricaviamo un Jacobiano valido per l'uso nel risolvere un vincolo di distanza. Possiamo iniziare dal passaggio 1 con \ eqref eq2. Ecco alcuni dettagli per la fase 2:

\ Begin equation
\ dot C = (P_2 - P_1) (\ dot P _2 - \ dot P _1)
\ Label eq9
\ End equation

\ Begin equation
\ dot C = (P_2 - P_1) ((v_2 + ω_2 \ volte r_2) - (v_1 + ω_1 \ volte r_1))
\ Label EQ10
\ End equation

\ (r_1 \) e \ (r_2 \) sono vettori dal centro di massa al punto di ancoraggio, rispettivamente per i corpi 1 e 2.

Il prossimo passo è quello di isolare i termini di velocità. Per fare ciò, utilizzeremo l'identità scalare del triplo prodotto:

\ Begin equation
(P_2 - P_1) = d
\ Label eq11
\ End equation

\ Begin equation
\ dot C = (d \ cdot v_2 + d \ cdot ω_2 \ times r_2) - (d \ cdot v_1 + d \ cdot ω_1 \ times r_1)
\ Label eq12
\ End equation

\ Begin equation
\ dot C = (d \ cdot v_2 + ω_2 \ cdot r_2 \ times d) - (d \ cdot v_1 + ω_1 \ cdot r_1 \ times d)
\ Label eq13
\ End equation

L'ultimo passo è identificare il Jacobiano mediante ispezione. Per fare ciò, tutti i coefficienti di tutti i termini di velocità (\ (V \) e \ (ω \)) saranno usati come elementi di Jacobian. Perciò:

\ Begin equation
J = \ begin bmatrix -d & -r_1 \ times d & d & r_2 \ times d \ end bmatrix
\ Label eq14
\ End equation

Altri Jacobians

Vincolo di contatto (vincolo di compenetrazione), dove \ (n \) è il contatto normale:

\ Begin equation
J = \ begin bmatrix -n & -r_1 \ times n & n & r_2 \ times n \ end bmatrix
\ Label eq15
\ End equation

Vincolo di attrito (attivo durante la penetrazione), dove \ (t \) è un asse di attrito (il 2D ha un asse, il 3D ne ha due):

\ Begin equation
J = \ begin bmatrix -t & -r_1 \ times t & t & r_2 \ times t \ end bmatrix
\ Label eq16
\ End equation


Risoluzione dei vincoli

Ora che abbiamo una comprensione di ciò che è un vincolo, possiamo parlare di come risolverli. Come affermato in precedenza, una volta che un Jacobiano è derivato dalla mano, dobbiamo solo risolvere per \ (\ lambda \). Risolvere un singolo vincolo in isolamento è facile, ma risolvere molti vincoli contemporaneamente è difficile e molto inefficiente (dal punto di vista computazionale). Questo pone un problema, poiché i giochi e le simulazioni probabilmente vorranno avere molti vincoli attivi contemporaneamente.

Un metodo alternativo per risolvere tutti i vincoli simultaneamente (risolvere globalmente) sarebbe quello di risolvere i vincoli in modo iterativo. Risolvendo per approssimazioni della soluzione e alimentando soluzioni precedenti alle equazioni, possiamo convergere sulla soluzione.

Uno di questi solutori iterativi è noto come Sequential Impulses, come soprannominato da Erin Catto. Sequential Impulses è molto simile a Projected Gauss Seidel. L'idea è di risolvere tutti i vincoli, uno alla volta, più volte. Le soluzioni si annulleranno a vicenda, ma su molte iterazioni ogni singolo vincolo convergerà e una soluzione globale può essere raggiunta. Questo è buono! I risolutori iterativi sono veloci.

Una volta raggiunta una soluzione, è possibile applicare un impulso a entrambi i corpi nel vincolo per far rispettare il vincolo.

Per risolvere un singolo vincolo, possiamo usare la seguente equazione:

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

\ (M ^ - 1 \) è la massa del vincolo; \ (b \) è il pregiudizio (più su questo più tardi).

Questa è una matrice contenente la massa inversa e l'inerzia inversa di entrambi i corpi rigidi nel vincolo. Quanto segue è la massa del vincolo; nota che \ (m ^ - 1 \) è la massa inversa di un corpo, mentre \ (I ^ - 1 \) è l'inerzia inversa di un corpo:

\ begin equation M ^ - 1 =
\ Begin bmatrix
m_1 ^ - 1 e 0 & 0 & 0 \\
0 e I_1 ^ - 1 e 0 e 0 \\
0 & 0 e m_2 ^ - 1 e 0 \\
0 & 0 & 0 & I_2 ​​^ - 1
\ End bmatrix
\ Label eq18
\ End equation

Anche se \ (M ^ - 1 \) è teoricamente una matrice, per favore non modellarla come tale (la maggior parte è zero!). Invece, sii intelligente su che tipo di calcoli fai.

\ (JM ^ - 1 J ^ T \) è noto come vincolo di massa. Questo termine è calcolato una volta e utilizzato per risolvere \ (\ lambda \). Lo calcoliamo per un sistema come questo:

\ Begin equation
JM ^ - 1 J ^ T = (l_1 \ cdot l_1) * m_1 ^ - 1 + (l_2 \ cdot l_2) * m_2 ^ - 1 + a_1 * (I_1 ^ - 1 a_1) + a_2 * (I_2 ^ - 1 a_2)
\ Label eq19
\ End equation

Si noti che è necessario invertire \ eqref eq19 per calcolare \ eqref eq17.

Le informazioni di cui sopra sono tutto ciò che è necessario per risolvere un vincolo! Una forza nello spazio cartesiano può essere risolta e utilizzata per aggiornare la velocità di un oggetto, al fine di imporre un vincolo. Si prega di richiamare \ eqref eq5:

\ Begin equation
F = \ lambda J ^ T \\
V_ final = V_ iniziale + m ^ - 1 * F \\
∴ \\
\ Begin bmatrix
v_1 \\
ω_1 \\
v_2 \\
ω_2 \\
\ end bmatrix + = \ begin bmatrix
m_1 ^ - 1 e 0 & 0 & 0 \\
0 e I_1 ^ - 1 e 0 e 0 \\
0 & 0 e m_2 ^ - 1 e 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 equation


Deriva dei vincoli

A causa della linearizzazione delle equazioni di posizione non lineari, alcune informazioni vengono perse. Ciò si traduce in soluzioni che non soddisfano completamente l'equazione di posizione originale, ma soddisfano le equazioni di velocità. Questo errore è noto come deriva dei vincoli. Si può pensare a questo errore come il risultato di un'approssimazione della linea tangente.

Ci sono alcuni modi diversi per risolvere tali errori, ognuno dei quali approssima l'errore e applica qualche forma di correzione. Il più semplice è noto come Baumgarte.

Baumgarte è una piccola aggiunta di energia nello spazio dei vincoli e tiene conto del termine \ (b \) nelle equazioni precedenti. Per tenere conto del bias, ecco una versione modificata di \ eqref eq4:

\ Begin equation
\ dot C = JV + b = 0
\ Label eq21
\ End equation

Per calcolare un termine Baumgarte e applicarlo come bias, dobbiamo esaminare l'equazione del vincolo originale e identificare un metodo adatto per calcolare l'errore. Baumgarte è nella forma:

\ Begin equation
JV = - \ beta C
\ Label eq22
\ End equation

\ (\ beta \) (termine Baumgarte) è un fattore sintonizzabile, dipendente dalla simulazione e senza unità. Di solito è tra 0.1 e 0.3.

Per calcolare il termine di bias, osserviamo l'equazione per il vincolo di non penetrazione \ eqref eq15 prima di derivare rispetto al tempo, dove \ (n \) è il contatto normale:

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

L'equazione di cui sopra sta dicendo che l'errore scalare di \ (C \) è la profondità di inter-penetrazione tra due corpi rigidi.


Conclusione

Grazie a Erin Catto e al suo articolo sulla risoluzione dei vincoli, abbiamo un modo per risolvere i vincoli di posizione in termini di velocità. Molti comportamenti interessanti possono derivare da vincoli e, si spera, l'articolo sarà utile a molte persone. Come sempre, non esitate a fare domande o dare commenti di seguito.

Si prega di vedere Box2D Lite per una risorsa sulla risoluzione di vari tipi di vincoli 2D, nonché informazioni dettagliate su molte specifiche di implementazione non trattate qui.