Numerisch stabile Methode zur Berechnung von Winkeln zwischen Vektoren

13

Bei Anwendung der klassischen Formel für den Winkel zwischen zwei Vektoren:

α=arccosv1v2v1v2

man stellt fest, dass bei sehr kleinen / spitzen Winkeln ein Präzisionsverlust auftritt und das Ergebnis nicht genau ist. Wie in dieser Antwort zum Stapelüberlauf erläutert , besteht eine Lösung darin, stattdessen den Arkustangens zu verwenden:

α=arctan2(v1×v2,v1v2)

Und das ergibt in der Tat bessere Ergebnisse. Ich frage mich jedoch, ob dies bei Winkeln nahe \ pi / 2 zu schlechten Ergebnissen führen würde π/2. Ist es der fall Wenn ja, gibt es eine Formel, mit der Winkel genau berechnet werden können, ohne eine Toleranz innerhalb eines ifZweigs zu prüfen ?

Astrojuanlu
quelle
1
Dies hängt von der Implementierung der inversen Tangensfunktion mit zwei Parametern ab. Die langsamen, stabilen Versionen wechseln bedingt zwischen der Arbeit mit x / y und y / x, um die Präzision beizubehalten, während die schnellen nur den richtigen Quadranten einhalten und daher nicht genauer sind als die Version mit einem Parameter.
Origimbo
Sie sollten "Genauigkeitsverlust" definieren: Angenommen, die richtige Antwort lautet und Sie erhalten stattdessenα + Δ Δ α Δ παα+Δ . Benötigen Sie oder ist ausreichend? ΔαΔπ
Stefano M
In diesem Fall war die richtige Antwort und ich habe , beide . α 10 81αα1081
Astrojuanlu

Antworten:

18

( Ich habe diesen Ansatz bereits getestet und erinnere mich, dass er korrekt funktioniert hat, aber ich habe ihn nicht speziell für diese Frage getestet. )

Soweit ich das beurteilen kann, sind beide und können unter einer katastrophalen Löschung leiden, wenn sie nahezu parallel / senkrecht sind - atan2 kann keine gute Genauigkeit liefern, wenn einer der Eingänge ausgeschaltet ist.v 1v 2v1×v2v1v2

Beginnen Sie mit einer Neuformulierung des Problems, indem Sie den Winkel eines Dreiecks mit den Seitenlängen,und(Diese werden alle in Gleitkomma-Arithmetik genau berechnet.) Es gibt eine bekannte Variante der Heronschen Formel aufgrund von Kahan ( Fehlberechnungsfläche und Winkel eines nadelförmigen Dreiecks ), mit der Sie die Fläche und den Winkel (zwischen und ) eines Dreiecks berechnen können, die durch seine Seitenlängen angegeben sind. und dies numerisch stabil. Da die Reduktion auf dieses Unterproblem ebenfalls genau ist, sollte dieser Ansatz für beliebige Eingaben funktionieren.b = | v 2 | c = | v 1 - v 2 |ein=|v1|b=|v2|c=|v1v2|bab

Zitiert aus diesem Artikel (siehe S.3), unter der Annahme, dass , Alle Klammern hier sind sorgfältig gesetzt, und sie sind wichtig; Wenn Sie feststellen, dass Sie die Quadratwurzel einer negativen Zahl ziehen, sind die eingegebenen Seitenlängen nicht die Seitenlängen eines Dreiecks.μ = { c - ( a - b ) , wenn  b c 0 , b - ( a - c ) , wenn  c > b 0 , ungültig Dreieck , sonst a n g l e = 2 arctan ( ab

μ={c(ab),if bc0,b(ac),if c>b0,invalid triangle,otherwise
angle=2arctan(((ab)+c)μ(a+(b+c))((ac)+b))

Wie dies funktioniert, wird in Kahans Aufsatz erklärt, einschließlich Beispielen für Werte, für die andere Formeln versagen. Ihre erste Formel für lautet auf Seite 4.αC

Der Hauptgrund, warum ich Kahans Herons Formel vorschlage, ist, dass sie ein sehr nettes Primitiv darstellt - viele potenziell knifflige Fragen der planaren Geometrie können auf das Finden der Fläche / des Winkels eines beliebigen Dreiecks reduziert werden eine schöne stabile Formel dafür, und es ist nicht nötig, sich selbst etwas auszudenken.

Bearbeiten Nach Stefanos Kommentar habe ich eine Darstellung des relativen Fehlers für , ( Code ) erstellt. Die zwei Linien sind die relativen Fehler für und ;, wobei entlang der horizontalen Achse verläuft. Es scheint, dass es funktioniert. v1=(1,0)v2=(cosθ,sinθ)θ=ϵθ=π/2ϵϵBildbeschreibung hier eingeben

Kirill
quelle
Danke für den Link und die Antwort! Leider erscheint die zweite Formel, die ich geschrieben habe, nicht im Artikel. Andererseits kann diese Methode etwas komplex werden, da eine Projektion in 2D erforderlich ist.
Astrojuanlu
2
@astrojuanlu Hier gibt es keine Projektion auf 2d: Unabhängig von den beiden 3D-Vektoren definieren sie ein einziges (planares) Dreieck zwischen sich - Sie müssen nur die Seitenlängen kennen.
Kirill
Sie haben recht, mein Kommentar macht keinen Sinn. Ich habe in Koordinaten statt in Längen gedacht. Danke noch einmal!
Astrojuanlu
2
@astrojuanlu Eine weitere Sache , ich anmerken möchte: Es scheint , dass es einen formalen Beweis ist , dass der Bereich Formel in genau Wie Berechne die Fläche eines Dreiecks: a Formal Revisit , Sylvie Boldo , mit Flocq.
Kirill
cc<ϵMindest(ein,b)(v1-v2)
7

Die effiziente Antwort auf diese Frage findet sich nicht allzu überraschend in einer weiteren Anmerkung von Velvel Kahan :

α=2arctan(v1v1+v2v2,v1v1-v2v2)

arctan(x,y)(x,y)

(Ich habe eine Mathematica Demonstration von Kahan Formel hier .)

JM
quelle
arctan2
1
arctan(x,y)ATAN2(Y, X)