Hi, mangels FPU mache ich die Berechnung von Quadratwurzeln selbst: - die Zahl "m", aus der die Wurzel gezogen werden soll, wird um 4 Bit rechts geshiftet - das Ergebnis wird in einer Schleife mit sich selbst multipilizert und dann (per sukzessiver Approximation) abhängig davon, ob es größer oder kleiner als "m" ist immer so weit verändert, dass es sich "m" weiter annähhert. Nach etwa 10 Schleifendurchläufen breche ich dann ab und verwende das Ergebnis. Jetzt ließe sich die Genauigkeit noch deutlich erhöhhen (bzw. die Anzahl der maximalen Schleifendurchläufe verringern) wenn ich mit der initialen Ermittlung des Startwertes näher am Ergebnis dran wäre als meine Division durch 16. Wie könnte ich hier vorgehen um mit meinem Startwert genauer zu sein?
Was für Datentypen verwendest du? Was ist der (übliche) Wertebereich für die Eingabe? Wie führst du die sukzessive Approximation durch? Auf welchem Prozessor soll das schnell laufen?
Daniel schrieb: > Was für Datentypen verwendest du? 32 Bit signed > Was ist der (übliche) Wertebereich für die Eingabe? 10240..1024000000, eine Häufung für einen bestimmten, kleineren Bereich ist nicht vorhersagbar. > Wie führst du die sukzessive Approximation durch? Grob vereinfacht so: r=m>>4; // -> initialen Wert ermitteln for (i=0; i<10; i++) { inter=r*r; if (inter>m) r_max=r; else if (inter<m) r_min=r; else return r<<5; r=(r_max+r_min)>>1; } return r<<5; > Auf welchem Prozessor soll das schnell laufen? ARM Cortex M3 - am liebsten aber Portabel in reinem C, ohne Assembler
Beim Cortex M3 hast du natürlich einen Vorteil durch CLZ (sicher auch als Intrinsic in deinem C/C++ Compiler verfügbar) Dann kannst du den Startwert besser annähern als mit /16 in dem du dir klar machst, dass bei der Multiplikation zweier Zahlen sich deren effektive Länge ca verdoppelt. Die effektive Länge von m ist 32-CLZ(m) die Länge von Wurzel aus m ist daher nicht kleiner als LowerGaussian((32-CLZ(m)) / 2) und nicht größer als LowerGaussian((32-CLZ(m)) / 2)+1 Du kannst also deinen Wertebereich auf 1 << (((32-CLZ(m)) / 2)-1) bis 1 << ((32-CLZ(m)) / 2) einschränken. Starte in der Mitte, dann hast du durchschnittlich den den "kürzeszten" Weg. 3 << (((32-CLZ(m)) / 2)-2) (hoffe ich hab keine allzugroßen Typos im Schnelltipp)
CLZ: _C_ount _l_eading _z_eros Aber was genau ist LowerGaussian() und aus welcher Bibliothek kommt die Funktion?
Lower Gaussian ist die mathematische Bezeichnung für die Rundung auf die nächstgelegene Ganzzahl mit kleinerem Absolutwert.
Man könnte auch einfach m>>((32-CLZ(m))>>1) als erste Näherung nehmen.
Ist bescheiden. zwar beliebig wählbar genau (inklusive einer garantierten Genauigkeit nach x Schritten) und nciht auf die Quadratwurzel beschränkt, aber hier zu aufwendig. Riesiger Minuspunkt: die Division.
öööhm... jemand itneressiert an einer exakten Fixkomma-Sqrt in 176 Takten auf nem AVR32..? Ist leider TemplateMagie in C++ -> aber um das Prinzip zu verstehen.. sollte in C konvertierbar sein. Der Algo benötigt genau 1 Division.
...hab mir damals ein PDF darüber gemacht, weil ich es für wirklich sehr schnell hielt^^ Code geb ich ungern raus... zumal er den Anforderungen nicht entspricht. (->kein C) Außerdem wird man nur bekloppt vollgemeckter hier... blabla. kein bock auf Codekritik.
DerAlbi schrieb: > Außerdem wird man nur bekloppt vollgemeckter hier... blabla. kein bock > auf Codekritik. Dann lass es doch. Oder geh in ein Forum, wo man sagt was man tolles hat, es aber nicht zeigt.
-.- Typ - Datentyp (un)signed int/short Fracbits: Anzahl der Kommastellen-bits Multiplikation mit floats -> bitte händisch durch fixkommarechung ersetzten (macht sonst meine lib) Fixkommat-Typ-Rechungen und Initialisierungen bitte händisch umändern. mac -> multiply accumulate Mac(a,b,c) = a*b+c div64u... nuja.. damit muss man selbst klarkommen. ich hab eine Divisionsroutine geschrieben die mit einem Hardeware-32Bit-Div eine gute näherung des Ergebnisses bringt. Was der Compiler bei einem echten 64Div einsetzt... dauert so lange wie die wurzel selbst^^
DerAlbi schrieb: > ...hab mir damals ein PDF darüber gemacht Sorry, du hast in entscheidenden Punkt einen Fehler: #1 Heronformel 2mal ineinander eingesetzt. # bn+1 = 1/2 * (bn + x/bn) = (x^2 + bn) / (2bn) //Gleichnamig beim Umformen sind dir x und bn durcheinander geraten. Es müsste heissen = (bn^2 + x) / 2bn Multiplikation von bn mit bn/bn(=1), und Addition beider Brüche. Das verhindet deine Vereinfachung auf einen Bruch ohne Abhängigkeit vom vorhergehenden Iterationsschritt.
Das is echt lachhaft :-D Irgendwie haste recht.... :-D Kann jemand erklären warums trotzdem funktioniert? :-)
Das tut es nicht. Beispiel bei dem das veränderte Verfahren fehlverhält: x:=16, Startwert: bn:=8 -> bn+1 = 16.5 -> bn+2 = 8.25 -> bn+3 = 16.02 -> bn+4 = 8.49 ... Es konvergiert nicht gegen 4
Ich sterb gleich. Was bitte ist denn da jetzt los. Ich hab damit nen 3D Renderer gebaut - und das Bild stimmt^^ DAS KANN doch gar nicht falschn sein... HÄÄÄÄÄÄÄÄÄÄÄÄÄÄ ..ich arbeite dran^^
Fehler gefunden :-) die Formel mit dem gemeinsamen Nenner ist falsch. Also die letzte Zeile. das was rauskommen soll lautet: ( x^2 + 6*bn^2*x + bn^4 ) / (4*bn*x + 4*bn^3) BEISPIEL: Wurzel aus 14 ld(14) = 3. ==> bn = 2^1.5 (14^2 + 48*14 + 64) / (4 * 2^1.5 * 14 + 4*2^4.5) = 3.744 Realer wert: 3.742 Ich weiß echt nicht, warum es funktioniert :-)
Also ich kanns immernoch nicht nachvollziehen :-/ Wie bin ich denn da drauf gekommen... die Heron-Formel ist ursprünglich korrekt, und wie sich herausgestellt hat... http://en.wikipedia.org/wiki/Methods_of_computing_square_roots#Bakhshali_approximation ist meine Endformel auch korrekt. Wie kommt man da auf mathematisch nachvollziehbaren Weg hin? :-D
Bitte melde dich an um einen Beitrag zu schreiben. Anmeldung ist kostenlos und dauert nur eine Minute.
Bestehender Account
Schon ein Account bei Google/GoogleMail? Keine Anmeldung erforderlich!
Mit Google-Account einloggen
Mit Google-Account einloggen
Noch kein Account? Hier anmelden.