Forum: Mikrocontroller und Digitale Elektronik Schneller sqrt()-Ersatz ohne FPU


von Floro (Gast)


Lesenswert?

Hi,

ich versuche gerade, Code auf einen ARM-Mikrocontroller ohne FPU zu 
portieren, der im Original vom PC kommt. Sämtliche 
Fließkommaberechnungen kann ich ausreichend genau ersetzen, in dem ich 
alle Werte als int behandle und einfach um den Faktor 1000 größer mache.

Einzig verbleibendes Problem: die Funktion sqrt(). Gibt es dafür 
irgendwo einen brauchbar schnellen Ersatzcode? Der mögliche Wertebereich 
liegt bei 30..1966020, so dass vorberechnete Tabellen leider wegfallen.

Irgend welche Ideen?

Danke!

von Karl H. (kbuchegg)


Lesenswert?


: Bearbeitet durch User
von Timm T. (Gast)


Angehängte Dateien:

Lesenswert?

Der beste Algo, den ich mal für Int-Sqrt auf dem µC gefunden habe, ist 
hier http://xlinux.nist.gov/dads//HTML/squareRoot.html beschrieben. Ist 
relativ schnell und gut erweiterbar, also 32bit sind kein Problem, 
bleiben ja eh nur 16bit nach der Ausführung übrig.

Für 32bit heisst schnell etwa 1000 Zyklen. Drunter wird es Dir auch ein 
anderer Algo ohne Tabelle nicht machen.

Ein Beispiel für PIC gibt es in AN040 von Microchip. Das kann man auf 
einen ordentlichen Controller portieren.

von Karl H. (kbuchegg)


Lesenswert?

Floro schrieb:

> Fließkommaberechnungen kann ich ausreichend genau ersetzen, in dem ich
> alle Werte als int behandle und einfach um den Faktor 1000 größer mache.

Das könnte ein Problem werden, das du berücksichtigen musst.
Denn die sqrt von 1000 muss ja logischerweise wieder 1000 ergeben.

von Floro (Gast)


Lesenswert?

Karl Heinz schrieb:
> Das könnte ein Problem werden, das du berücksichtigen musst.
> Denn die sqrt von 1000 muss ja logischerweise wieder 1000 ergeben.

Ja, das ist in dem zu erwartenden Wertebereich bereits enthalten.

von Martin L. (maveric00)


Lesenswert?

Hallo,

je nachdem, wieviel Flash Du noch frei hast, könnte eine Tabelle 
durchaus möglich sein, sie würde rund 5,5 kB benötigen.

Die Grundidee wäre dann, bei y=sqrt(x) keine Tabelle für alle "x" mit 
Verweis auf "y" anzulegen, sondern eine Tabelle mit den Werten von y^2, 
wobei y der Index wäre und dann mit einem Suchalgorithmus die beiden 
Indizes zu finden, für die y1^2<x<y2^2 gilt. Danach kann man 
entsprechend auf y1 ab- oder auf y2 aufrunden. Sollte ziemlich schnell 
sein (jedenfalls deutlich weniger als 1000 Zyklen brauchen).

Schöne Grüße,
Martin

von Rolf Magnus (Gast)


Lesenswert?

Floro schrieb:
> Ja, das ist in dem zu erwartenden Wertebereich bereits enthalten.

Das Näherungsverfahren muß den Faktor aber auch berücksichtigen. Wie 
Karl Heinz schon schreibt:

Karl Heinz schrieb:
> Denn die sqrt von 1000 muss ja logischerweise wieder 1000 ergeben.

Ohne Berücksichtigung des Faktors 1000 käme ca.31 raus statt 1000. Damit 
wäre die Wurzel von 1 nicht 1, sondern 0,031.

von Karl H. (kbuchegg)


Lesenswert?

Floro schrieb:
> Karl Heinz schrieb:
>> Das könnte ein Problem werden, das du berücksichtigen musst.
>> Denn die sqrt von 1000 muss ja logischerweise wieder 1000 ergeben.
>
> Ja, das ist in dem zu erwartenden Wertebereich bereits enthalten.

SChon klar.
Aber waruf ich raus will ist, dass dir ein normaler Wurzelalgorithmus 
als sqrt(1000) das Ergebnis 31 liefern wird (oder 32 wenn er rundet).
Das kannst du aber nicht brauchen. sqrt(1000) muss 1000 als Ergebnis 
haben.


Ich würde mein Heil zuerst mal wahrscheinlich in der Digit-Methode 
suchen. Deine int sind 32 Bit groß? dann ist klar, dass das Ergebnis 
nicht größer als 16 Bit sein kann.
Nacheinander dann, beginnend vom MSB eines 16 Bit Wertes (den ich mal x 
nenne), testweise ein 1 Bit setzen, die Kontrollmultiplikation machen 
und nachsehen ob x*x größer als die Zahl ist, von der die Wurzel gezogen 
werden soll. Wenn größer, dann das Bit wieder auf 0 setzen und das 
nächstkleinere Bit betrachten. Durch deine Spezialmultiplikation, die 
die besonderheiten der Fixed Point Arithmetik schon berücksichtigt, 
sollte da dann auch das richtige rauskommen

pseudo
1
int fix_sqrt( int a )
2
{
3
  int result;
4
  int i;
5
6
  for( i = 0; i < 16; ++i )
7
  {
8
    result |= ( 1 << ( 15 - i ) );
9
    if( fix_mult( result, result ) > a )
10
      result &= ~( 1 << (15 - i ) );
11
  }
12
13
  return result;
14
}

ungetester Code, aber so ungefähr müsste das klappen.


PS: die 1000 waren keine gute Idee. DU halst damit deinem µC mehr Arbeit 
auf als notwendig. 10-er Potenzen sind zwar für uns Menschen intuitiv, 
aber das was für uns 10-er Potenzen an Vereinfachungen bringen, das sind 
für einen Computer 2-er Potenzen. 1024 wäre programmtechnisch besser 
gewesen, auch wenn die Zahlen dann im Debugger anders aussehen.

: Bearbeitet durch User
von Lothar M. (Firma: Titel) (lkmiller) (Moderator) Benutzerseite


Lesenswert?

Floro schrieb:
> Gibt es dafür irgendwo einen brauchbar schnellen Ersatzcode?
Siehe dort (ganz unten):
http://www.lothar-miller.de/s9y/archives/73-Wurzel-in-VHDL.html

von ARMdran (Gast)


Lesenswert?

Blöde Frage, aber hast Du mal versucht, den Code ohne Modifikationen 
hinsichtlich fixed point zu portieren?

Man wundert sich manchmal, was so ein kleiner ARM (auch ohne FPU) mit 
entsprechender float-lib so leistet. Ich würde auf jeden Fall bei floats 
bleiben, solange es irgendwie geht.

von Markus M. (Firma: EleLa - www.elela.de) (mmvisual)


Lesenswert?

Welcher ARM uC?

von Marc P. (marcvonwindscooting)


Lesenswert?

Floro schrieb:
> Sämtliche
> Fließkommaberechnungen kann ich ausreichend genau ersetzen, in dem ich
> alle Werte als int behandle und einfach um den Faktor 1000 größer mache.

Ich w"urde die Werte alle um 1024 gr"osser machen.

@ARMdran
Ich w"urde so fr"uh wie m"oglich floats vermeiden, denn eine Addition in 
einem Takt und KEIN float-lib aufruf ist schon ein gewaltiger 
Unterschied.
Mich wundert manchmal, was so ein kleiner ARM (LPC176x) ohne FPU und 
ohne floats so leisten kann. Zum Beispiel in einem 3P-Frequenzumrichter 
mit bis zu 70000 samples / Sekunde 3 Kinderpopo-glatte Sinuskurven 
(768level, 10bit-Zeit, DDS) mit Amplitude, 
Motorbeschleunigungskennlinie, Spannungskennlinie in einer ISR zu machen 
und DOCH NOCH ANDERE Sachen tun...

An Festkomma muss man sich echt gew"ohnen, aber daf"ur geht der Proz 
dann ab dass man beim Messen der Geschwindigkeit fast aus den Latschen 
kippt.

von ARMdran (Gast)


Lesenswert?

Klar ist fixed point schneller, aber ich seh das so:

Man designed ein Gerät ja nach den Anforderungen. Und wenn man die 
Anforderungen mit floats heutzutage nicht erreichen kann, hat man mit 
großer Wahrscheinlichkeit den falschen Prozessor selektiert.
Fixed points zu verwenden nur weil es schneller ist (obwohl man die 
Systemanforderungen mit floats auch erreichen könnte) macht keinen Sinn.

von Marc P. (marcvonwindscooting)


Lesenswert?

ARMdran schrieb:
> Fixed points zu verwenden nur weil es schneller ist (obwohl man die
> Systemanforderungen mit floats auch erreichen könnte) macht keinen Sinn.

Doch. Macht es.
Es macht den Sinn der geringeren Materialkosten, des geringern 
Strombedarfs, des geringern Speicherbedarfs ( = schon wieder Strombedarf 
und Kosten).

Es macht den Sinn, dass manche Integer-Arithmetik genauer ist als 
floats. Schau die mal die mickrige Mantisse von 'ne 4-byte float an. 
Okay dann halt double und der n"achste Prozessor bitte...

Wir sind heutzutage st"andig versucht, alle Prozessorleistung mit 
Schnickschnack abzufackeln, nur damit nix "ubrig bleibt um mal die 
eigentliche Aufgabe BESSER zu l"osen als das fr"uher ging.

von Robin (Gast)


Lesenswert?

Rolf Magnus schrieb:
> Ohne Berücksichtigung des Faktors 1000 käme ca.31 raus statt 1000. Damit
> wäre die Wurzel von 1 nicht 1, sondern 0,031.

Wenn man das ergebniss dann mal Wurzel(1000) rechnet kommen am ende auch 
wieder ca. 1000 raus. Das problem sind halt die Kommastellen, die 
verloren gehen. Besser wäre wohl ein anderer Faktor, der auch eine 
Ganzzahlige wurzel hat, dann ist der Fehler am ende kleiner.

von ARMdran (Gast)


Lesenswert?

Marc P. schrieb:
> ARMdran schrieb:
>> Fixed points zu verwenden nur weil es schneller ist (obwohl man die
>> Systemanforderungen mit floats auch erreichen könnte) macht keinen Sinn.
>
> Doch. Macht es.
> Es macht den Sinn der geringeren Materialkosten, des geringern
> Strombedarfs, des geringern Speicherbedarfs ( = schon wieder Strombedarf
> und Kosten).

Entwickelst du in Stückzahlen > 1000? Ok, dann mag das Sinn machen. Der 
Speicherbedarf geht kaum in den Strombedarf ein.


> Es macht den Sinn, dass manche Integer-Arithmetik genauer ist als
> floats. Schau die mal die mickrige Mantisse von 'ne 4-byte float an.
> Okay dann halt double und der n"achste Prozessor bitte...

Das mag für bestimmte Fälle gelten. Dafür sind die Entwicklungskosten 
für die fixed-point-Entwicklung i.d.R höher. Und nicht jeder 
Codegenerator beherrscht fixed points sinnvoll.

>
> Wir sind heutzutage st"andig versucht, alle Prozessorleistung mit
> Schnickschnack abzufackeln, nur damit nix "ubrig bleibt um mal die
> eigentliche Aufgabe BESSER zu l"osen als das fr"uher ging.

Eine Aufgabe ist dann gelöst, wenn die Requirements erfüllt sind.

von asd (Gast)


Lesenswert?

Z.B. 1024 und 32. Damit wird das Multiplizieren und Dividieren zu 
einfachen Shifts.

von Marc P. (marcvonwindscooting)


Lesenswert?

Wenn Du '/' hast, dann kannst Du auch mal ganz grob sowas machen:

Eingabe: x

w = x/2;
for (int i=0; i<Geduld; i++) w = (w + x/w) / 2

Ausgabe w.


Die alte Rekursionsformel zieht auch ganz sch"on schnell an....

von W.S. (Gast)


Lesenswert?

Floro schrieb:
> Sämtliche
> Fließkommaberechnungen kann ich ausreichend genau ersetzen, in dem ich
> alle Werte als int behandle und einfach um den Faktor 1000 größer mache.

Wie bitte???
Das funktioniert nur mit Addition und Subtraktion. Bei Multiplikation 
und Division wirst du dein blaues Wunder erleben. Wenn schon, dann 1024 
und Shifts nach allen entsprechenden Operationen.

Ansonsten lies dir mal was zum Thema Pseudodivision und 
Pseudomultiplikation durch. Klug eingesetzt kann das ein netter Weg für 
transzendente Funktionen und Ähnliches sein.

W.S.

von Lothar M. (Firma: Titel) (lkmiller) (Moderator) Benutzerseite


Lesenswert?

Karl Heinz schrieb:
> 
http://www.drdobbs.com/embedded-systems/square-root-part-2-or-1414-if-you-prefer/229500627
Dankesehr, jetzt weiß ich endlich wieder, wo ich das her hatte...  ;-)

Marc P. schrieb:
> Mich wundert manchmal, was so ein kleiner ARM (LPC176x) ohne FPU und
> ohne floats so leisten kann. Zum Beispiel ...
Ja, da fragt man sich, wo die PC-Programmierer die vielen Cores und die 
Gigahertzen hintun  :-/

von Marc P. (marcvonwindscooting)


Lesenswert?

Marc P. schrieb:
> Die alte Rekursionsformel zieht auch ganz sch"on schnell an....

Geduld = 13 reicht in deinem Wertebereich.
Glaubt man gar nicht, dass die so schnell konvergiert, oder?

von Marc P. (marcvonwindscooting)


Lesenswert?

@Floro: vernachl"assig mal deinen Fredd nich so! ;-)

Was macht denn deine Festpunkt-Umsetzung? Welche Wurzel?

Marc P. schrieb:
> Ich w"urde die Werte alle um 1024 gr"osser machen.

ich m"ochte meine Aussage revidieren: Ich w"urde seit gestern die Werte 
alle um (den Faktor) 65536 (2^16) gr"osser machen.
Ich habe die vergangenen Tage u.a. meine Festpunkt-Ausgaberoutinen auf 
einen besseren Stand gebracht. Da ich auch den Cortex-M0 benutze muss 
ich in letzter Zeit wieder etwas mehr darauf achten, keine unn"otigen 
Divisionen zu machen. Und weil der M0 auch keinen 'long multiply' 
(SMULL, UMULL) hat, kann er nur 16x16 => 32 bit sehr schnell 
multiplizieren. So ergibt sich letztendlich, dass Festkommazahlen mit 
Komma zwischen bit15 und bit16 ab besten dezimal auszugeben sind - 4..5 
MAL so schnell wie andere: 25 700 Konvertierungen/s @ 12MHz Takt, um 
mal einen absoluten Wert in den Raum zu stellen. Das sind 467 Takte pro 
Konvertierung (5vor dem Komma, 4nach dem Komma, inclusive padding) mit 
der Luxus-Routine.

von Timm T. (Gast)


Lesenswert?

Marc P. schrieb:
> mal einen absoluten Wert in den Raum zu stellen. Das sind 467 Takte pro
> Konvertierung (5vor dem Komma, 4nach dem Komma, inclusive padding) mit
> der Luxus-Routine.

Du wandelst eine 32bit-Zahl in eine 5+4stellige Dezimalzahl um? Und das 
soll 500 Takte dauern?

Wie machst Du das?

von Arne (Gast)


Angehängte Dateien:

Lesenswert?

s.Anhang

von Marc P. (marcvonwindscooting)


Lesenswert?

Marc P. schrieb:
> w = x/2;
> for (int i=0; i<Geduld; i++) w = (w + x/w) / 2

Ist das Gleiche, nur mit einem anderen Startwert und weniger schwer 
verst"andlichen Zuweisungen und einem sehr kleinen Geduldsfaktor ;-)

Timm Thaler schrieb:
> Du wandelst eine 32bit-Zahl in eine 5+4stellige Dezimalzahl um? Und das
> soll 500 Takte dauern?
>
> Wie machst Du das?

Du verunsicherst mich. Entweder bin ich so grottenschlecht 
(50Takte/digit) und hab's bloss bisher nicht kapiert, oder Du findest es 
tats"achlich schnell!?

Lad dir mal bei:
http://www.windscooting.com/softy/software.html

das c-linux-r89.tar.bz2 runter und schau nach in:
c-linux/lib/c-any/{printFixedPoint.h,printFixedPoint-sliced.c,uint16Div. 
h}

Nimm an, dass USE_MX16DIV10 '#defined' ist, wenn Du f"ur einen Cortex-M0 
compilierst.

von Arc N. (arc)


Lesenswert?

Timm Thaler schrieb:
> Der beste Algo, den ich mal für Int-Sqrt auf dem µC gefunden habe, ist
> hier http://xlinux.nist.gov/dads//HTML/squareRoot.html beschrieben. Ist
> relativ schnell und gut erweiterbar, also 32bit sind kein Problem,
> bleiben ja eh nur 16bit nach der Ausführung übrig.
>
> Für 32bit heisst schnell etwa 1000 Zyklen. Drunter wird es Dir auch ein
> anderer Algo ohne Tabelle nicht machen.


16-Bit Wurzel aus einer 32-Bit Zahl...
"3 cycle/bit, 51 cycle total ARM code routine"
http://www.finesse.demon.co.uk/steven/sqrt.html

von Timm T. (Gast)


Lesenswert?

Arc Net schrieb:
> 16-Bit Wurzel aus einer 32-Bit Zahl...
> "3 cycle/bit, 51 cycle total ARM code routine"

Ich hab das damals auf dem ATmega8 gemacht. Da sind diese Befehle nicht 
verfügbar. Setzt man die Routine für den Mega8 um, dürfen auch ein paar 
mehr Zyklen rauskommen.

Ich die Rountine nochmal rausgekramt und getestet. Durchlauf für 32bit 
ist relativ unabhängig von der Größe der Zahl zwischen 250 und 280 
Zyklen. Die 1000 waren wohl etwas hoch geschätzt...

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.