Forum: Compiler & IDEs AVR: Floating Point => integer arithmetik


von Michael R. (Firma: Brainit GmbH) (fisa)


Lesenswert?

Hallo allerseits,

ich habe folgenden (rekursiven) Ausdruck, der mit float-Berechnugnen ca. 
600 Zyklen benötigt. Den würde ich nun gerne in integer-Arithmetik 
umwandeln, und könnte da ein paar Tipps gebrauchen:
1
p = p * (1 - m * p * p)
was weiss ich darüber:

p ist im Endeffekt ein uint16, und zwar der Maximalwert eines 
16Bit-Timers. ich gehe davon aus, dass der gesamte Bereich aber nicht 
ausgenützt wird (siehe aber weiter unten). Der Wert ist aber sicher 
größer als 1


Es gilt 0 < (1 - m*p*p) < 1
Das ist zwar gut zu wissen, aber für integer eher schlecht :-)

Daraus folgt aber auch dass 0 < m*p*p < 1

m wird berechnet, bleibt während der kritischen Schleifendurchläufe aber 
konstant

Der Ausdruck wurde (nicht von mir) so umgeformt, dass  er nur 
Multiplikationen enthält, aber keine Divisionen. Dass kann man aber 
ändern...

Der Ausdruck lässt sich auf verschiedene Arten umformen:

p *= (1-m*p*p) (hilft aber nicht viel)

p = p - m*p^3
p -= m*p^3

m kann durch den Kehrwert 1/n ersetzt werden:

p = p*(1-p*p/n)

aus obigem 0 < p*p/n < 1 folgt n > p*p

folgt daraus auch p^2 < n < p^3 ? da bin ich mir grad nicht sicher...

andere Gedanken:

um p^3 sicher in 32bit berechnen zu können, müsste p < Wurzel3(2^32) 
sein, also < 1625. Das wird eher eng... (gefühlsmäßig komme ich in den 
Bereich 10000..20000

Um die Berechnung möglichst genau auszuführen, könnte man vorher prüfen:

if p < 1625
   p*p*p/m (zuerst 2 mal multiplizieren dann dividieren)
else
   p*p/n*p (zuerst multiplizieren, dazwischen dividieren, dann nochmal 
multiplizieren

ist das sinnvoll? Muss man da befürchten dass einem der gcc da irgendwas 
dreht? (wie war das nochmal mit den sequence points?)

Edit: nein, das geht ja gar nciht, weil n > p*p, und damit wird p*p/n 
immer 0 ergeben :-(

Weiters: da gibts ja so eine komische Funktion, die einem bei 
integer-division in einem Rutsch auch den Rest liefert. Sollte man die 
verwenden? Wenn ja: wie lasse ich den Rest in die nächste Schleife 
einfließen?


Sonst noch Ideen zu meinem Problem?

Ach ja: Alternative mit Lookup-Table ist mir klar, aber wenns mit 
int-arithmetik einfacher ginge...

: Bearbeitet durch User
von Johann L. (gjlayde) Benutzerseite


Lesenswert?

Michael Reinelt schrieb:
> ich habe folgenden (rekursiven) Ausdruck, der mit float-Berechnugnen ca.
> 600 Zyklen benötigt. Den würde ich nun gerne in integer-Arithmetik
> umwandeln, und könnte da ein paar Tipps gebrauchen:

Wenn du mit 32-Bit int rechnest und Divisionen verwendest dürfte das 
eher langsamer werden als die 600 Ticks mit float.

> Daraus folgt aber auch dass 0 < m*p*p < 1

Woraus weiter folt daß m sehr klein ist; ohne Division kommst du also 
nicht aus.  Und wenn das Ergebnis ein int sein soll wird's in 0 < ... < 
1 recht eng ;-)

von Michael R. (Firma: Brainit GmbH) (fisa)


Lesenswert?

Johann L. schrieb:
> Wenn du mit 32-Bit int rechnest und Divisionen verwendest dürfte das
> eher langsamer werden als die 600 Ticks mit float.

Urgs... und ich hab mich noch gewundert, dass in der Doku zu dem 
Algorithmus stand "optimized for floating point, which will be faster 
than integer because it does not use divisions", hab das aber als Voodoo 
abgetan...

>> Daraus folgt aber auch dass 0 < m*p*p < 1
>
> Woraus weiter folt daß m sehr klein ist; ohne Division kommst du also
> nicht aus.  Und wenn das Ergebnis ein int sein soll wird's in 0 < ... <
> 1 recht eng ;-)

Ja, m ist sehr sehr klein.

Falls ich bei float bleibe: gibts sonst noch irgendwelche tricks, das zu 
beschleunigen? Fixed point wird ja auch nix bringen, oder?

von Falk B. (falk)


Lesenswert?

@ Michael Reinelt (fisa)

>Falls ich bei float bleibe: gibts sonst noch irgendwelche tricks, das zu
>beschleunigen? Fixed point wird ja auch nix bringen, oder?

Sage uns, was du INSGESAMT erreichen willst, siehe Netiquette.

von Michael R. (Firma: Brainit GmbH) (fisa)


Lesenswert?

Falk Brunner schrieb:
> Sage uns, was du INSGESAMT erreichen willst, siehe Netiquette.

Ich möchte den Algorithmus von Aryeh Eiderman implementieren: 
http://www.hwml.com/LeibRamp.htm Gleichung [20]

von Michael R. (Firma: Brainit GmbH) (fisa)


Lesenswert?

Johann L. schrieb:
> Wenn du mit 32-Bit int rechnest und Divisionen verwendest dürfte das
> eher langsamer werden als die 600 Ticks mit float.

Ich Ungläubiger hab das nun verifiziert: tatsächlich macht eine 
int32-Division den Code langsamer als das Multiplizieren mit floats... 
gut zu wissen!

von Detlef _. (detlef_a)


Lesenswert?

Hi,

das Problem ist die hohe Dynamik: Wenn p 16 Bit hat, hat p*p 32 Bit. Mit 
m*p*p~1 muß Du m mit 32 Bit skalieren. Ich denke, es ist allenfalls 
möglich den Ausdruck m*p*p mit zwei 16Bit*32Bit Multiplikationen zu 
rechnen.

Herr Eidermann macht aber auch Vorschläge zur besseren Berechenbarkeit 
der Gleichung. Ich würde mal ausprobieren, wie sich die Iterationen mit 
den Vereinfachungen verhalten, also mathematischer und nicht 
rechentechnisch agieren ;-)))))

math rulez
Cheers
Detlef

von Michael R. (Firma: Brainit GmbH) (fisa)


Lesenswert?

Detlef _a schrieb:
> Herr Eidermann macht aber auch Vorschläge zur besseren Berechenbarkeit
> der Gleichung. Ich würde mal ausprobieren, wie sich die Iterationen mit
> den Vereinfachungen verhalten, also mathematischer und nicht
> rechentechnisch agieren ;-)))))

ich bin ja schon bei seiner "einfachsten" Variante... oder hätte ich da 
was übersehen?

von Detlef _. (detlef_a)


Lesenswert?

Michael Reinelt schrieb:
> Detlef _a schrieb:
>> Herr Eidermann macht aber auch Vorschläge zur besseren Berechenbarkeit
>> der Gleichung. Ich würde mal ausprobieren, wie sich die Iterationen mit
>> den Vereinfachungen verhalten, also mathematischer und nicht
>> rechentechnisch agieren ;-)))))
>
> ich bin ja schon bei seiner "einfachsten" Variante... oder hätte ich da
> was übersehen?

Nee, Du hast Recht, ich hab übersehen, dass das die einfachste Form ist. 
Hab keine Idee.

Cheers
Detlef

von Bla (Gast)


Lesenswert?

Wenn p*p*p zu groß für uint32_t ist, dann geht vielleicht
1
uint32_t h = p * p;
2
h *= m;
3
h *= p;
4
p -= h;

Das produziert natürlich Fehler, je nach Wertebereich.

von Bla (Gast)


Lesenswert?

Wenn p*p*m < 1 geht das natürlich nicht.

m muss ja winzig sein! Woher kriegst du m? Kannst du m irgendwie 
speichern als Vielfaches von 1 / 2^32, also als ganze Zahl? Dann 
könntest du multiplizieren, shiften, multiplizieren, shiften...

von Robin (Gast)


Lesenswert?

Wie groß können deine Werte denn werden?

p liegt zwischen 1 und 65535
Aus 0 < m*p*p < 1 muss m im worst case (großes p) kleiner als 1/65536 
sein.

ich würde p = 16 bit, und m als 0.32 Festkomma behandeln.

p*p -> 32 bit
p*p*m -> bleibt 32 bit, da die oberen 32 Bit vernachlässigt werden 
können.

aus 1-p*p*m wird dann 65536-p*p*m und das Ergebniss ist immernoch ein 
0.32 Zahl.

Das ganze noch mal p rechnen und du hast im schlimmsten fall eine 16.32 
Zahl. Eine Division brauchst du da eigentlich nicht, musst nur mit 
ziemlich großen Zahlen rechnen.

von Quack (Gast)


Lesenswert?

Spannend wie viel man ueber die Sache diskutieren kann, ohne den 
Wertebereich von m zu kennen.

von Michael R. (Firma: Brainit GmbH) (fisa)


Lesenswert?

Bla schrieb:
> m muss ja winzig sein! Woher kriegst du m? Kannst du m irgendwie
> speichern als Vielfaches von 1 / 2^32, also als ganze Zahl? Dann
> könntest du multiplizieren, shiften, multiplizieren, shiften...

m ist winzig: m = a / F_CPU^2 mit a (Beschlunigung) 200.000 Steps/sec => 
m = 7.8E-10, * 2^32 = 3.35, und 3 ist als Ganzzahl wohl zu klein bzw. zu 
ungenau

Robin schrieb:
> ich würde p = 16 bit, und m als 0.32 Festkomma behandeln.
>
> p*p -> 32 bit
> p*p*m -> bleibt 32 bit, da die oberen 32 Bit vernachlässigt werden
> können.

Interessanter Ansatz! Nur kann ich nicht ganz folgen: wie mach ich das 
mit dem 0.32 Festkomma?

von Bla (Gast)


Lesenswert?

1
uint32_t h = 858; // 7.8E-10 * 2^40
2
h *= p;
3
h >>= 8;
4
h *= p;
5
h >>= 16;
6
h *= p;
7
h >>= 16;
8
p -= h;

von Robin (Gast)


Lesenswert?

Michael Reinelt schrieb:
>> p*p -> 32 bit
>> p*p*m -> bleibt 32 bit, da die oberen 32 Bit vernachlässigt werden
>> können.
>
> Interessanter Ansatz! Nur kann ich nicht ganz folgen: wie mach ich das
> mit dem 0.32 Festkomma?

Das ganze nennt sich Festkomma arithmetik, C kann das leider nicht von 
sich aus, d.h. du musst beim Berechnen selbst ans shiften denken und in 
welchem Zahlenformat sich deine Ergebnisse dann befinden. Aber du kannst 
recht einfach und schnell auf dem AVR mit Kommas rechnen.

von Johann L. (gjlayde) Benutzerseite


Lesenswert?

Robin schrieb:

> ich würde p = 16 bit, und m als 0.32 Festkomma behandeln.
>
> p*p -> 32 bit
> p*p*m -> bleibt 32 bit, da die oberen 32 Bit vernachlässigt werden
> können.

Was man vernachlässigen kann sind die unteren 32 Bit.  Die Rechnung 
müsste trotzdem in 64 Bit ausgeführt werden und dann nur die oberen 32 
Bits als Ergebnis genommen.  Oder eben ne eigene, spezielle 
Multiplikationsroutine schreiben, ist aber nicht so komfortabel zu 
implementieren wie eine 32*32 Multiplikation.

Michael Reinelt schrieb:
> ich bin ja schon bei seiner "einfachsten" Variante... oder hätte ich da
> was übersehen?

Wenn du bei float bleibst und m nur sehr selten zu berechnen ist, dann 
kann man eine Multiplikation sparen durch

Falls m nur die Werte -m0, 0 und m0 annimmt, kann man das einmal beim 
Programmstart vorberechnen.

Michael Reinelt schrieb:
> m ist winzig: m = a / F_CPU^2 mit a (Beschlunigung) 200.000 Steps/sec
> => m = 7.8E-10,

Der Motor macht 200000 Schritte pro Sekunde mit? Nicht wirklich, oder?

Wie oft wied der Berechnete Wert denn wirklich gebraucht? Evtl. ergibt 
dann Renormierung auf ein größeres Intervall besser handhabbare Werte?

von Michael R. (Firma: Brainit GmbH) (fisa)


Lesenswert?

Johann L. schrieb:
> Wenn du bei float bleibst und m nur sehr selten zu berechnen ist, dann
> kann man eine Multiplikation sparen durch
>
> Falls m nur die Werte -m0, 0 und m0 annimmt, kann man das einmal beim
> Programmstart vorberechnen.

Ähhh... klingt gut, aber wer war schnell q? Und du hast recht: m ist nur 
-R, 0 oder R

Nachdem ich zum Rasenmähen immer drei Bier brauche, steig ich da für 
heute aus :-) und versuchs morgen früh nochmal.


>
> Michael Reinelt schrieb:
>> m ist winzig: m = a / F_CPU^2 mit a (Beschlunigung) 200.000 Steps/sec
>> => m = 7.8E-10,
>
> Der Motor macht 200000 Schritte pro Sekunde mit? Nicht wirklich, oder?

Sorry, hab mich bei der Einheit vertan: Die Beschleunigung dürfte im 
bereich 200.000 Steps/sec/sec liegen, die Maximalgeschwindigkeit ist 
natürlich deutlich niedriger (aber immer noch recht hoch wegen 
Microstepping)


> Wie oft wied der Berechnete Wert denn wirklich gebraucht? Evtl. ergibt
> dann Renormierung auf ein größeres Intervall besser handhabbare Werte?

Tja, eigentlich nach jedem Step, zumindest in der Beschleunigungsphase.

von Michael R. (Firma: Brainit GmbH) (fisa)


Lesenswert?

Johann L. schrieb:
> Wenn du bei float bleibst und m nur sehr selten zu berechnen ist, dann
> kann man eine Multiplikation sparen durch

Zumindest die Formel habe ich verstanden, aber noch nicht wo ich mir die 
Multiplikation spare:

Lassen wir das Vorzeichen mal weg (Annahme: m immer positiv), dann kann 
ich mir sqrt(m) vorberechnen, nennen wir es n: n = sqrt(m)

Die neue Berechnung lautet dann: (n*q)^2 oder h=m*q, h*h.

ich seh aber immer noch zwei Multiplikationen?

von Johann L. (gjlayde) Benutzerseite


Lesenswert?

Oje, ich sollte echt keine Beiträge mehr schreiben, wenn ich schon total 
ausgenudelt bin...  Ich hatte im ursprüngichen Term 3 Multiplikationen 
gezählt obwohl es 3 Faktoren sind.

von Michael R. (Firma: Brainit GmbH) (fisa)


Lesenswert?

Johann L. schrieb:
> Ich hatte im ursprüngichen Term 3 Multiplikationen
> gezählt obwohl es 3 Faktoren sind.

Kein problem!

Aber jetzt musst du einem Österreicher nur noch erklären, was

Johann L. schrieb:
> ausgenudelt
ist ;-)

von Hans-jürgen H. (hjherbert) Benutzerseite


Angehängte Dateien:

Lesenswert?

Guten Abend

Ich sehe da keine Multiplikationen / Divisionen, außer mit Festwerten.
Und die kan man so wählen dass nur Bytes herumgehoben oder in 
schlechteren Falls durch 2er-Potenzen dividiert wird.

von Bla (Gast)


Lesenswert?

So hatte ich das oben irgendwo gepostet. Es sind Shift-Operationen statt 
Divisionen nötig, und in meinem Beispiel oben sind es nur 8 bzw. 16 Bit, 
d.h. der Compiler macht da auf einer 8-Bit-CPU gar keine Shift-Operation 
draus. Nimmt halt nur die Bytes für die nächsten Schritte, die er 
braucht. Das dürfe die Geschichte massiv ggü. float beschleunigen.

von Michael R. (Firma: Brainit GmbH) (fisa)


Lesenswert?

Bla schrieb:
> So hatte ich das oben irgendwo gepostet. Es sind Shift-Operationen statt
> Divisionen nötig, und in meinem Beispiel oben sind es nur 8 bzw. 16 Bit,
> d.h. der Compiler macht da auf einer 8-Bit-CPU gar keine Shift-Operation
> draus. Nimmt halt nur die Bytes für die nächsten Schritte, die er
> braucht. Das dürfe die Geschichte massiv ggü. float beschleunigen.

Yep, ich hab deinen Code dann eh aufgegriffen. Dein 2^40 hat mich zuerst 
verwirrt, aber nach etwas Nachdenken hab ichs kapiert.

von Bla (Gast)


Lesenswert?

16 + 16 + 8 = 40

Gern geschehen übrigens!

von Michael R. (Firma: Brainit GmbH) (fisa)


Lesenswert?

Bla schrieb:
> 16 + 16 + 8 = 40

Ja genau, den Schritt hab ich geistig erst nicht geschafft. 40 sieht so 
"unverdächtig" aus :-)

Schnell ist es tatsächlich: ich bin jetzt bei etwa 100 Takten (genau hab 
ichs nicht mehr im Kopf)


ah ja, nachträglich Danke, übrigens :-)

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
Noch kein Account? Hier anmelden.