Nächste: , Vorige: , Nach oben: Fourier-Transformationen   [Inhalt][Index]

23.2 Funktionen und Variablen für die schnelle Fourier-Transformation

Funktion: polartorect (r, t)

Transformiert komplexe Zahlen der Form r %e^(%i t) in die Standardform a + b %i. r ist der Betrag der komplexen Zahl und t die Phase. Die Argumente r und t sind eindimensionale Arrays derselben Größe. Die Größe der Arrays muss eine Potenz von 2 sein.

Die Werte der originalen Arrays werden durch den Realteil a = r cos(t) und den Imaginärteil b = r sin(t) ersetzt.

polartorect ist die inverse Funktion zu recttopolar. Das Kommando load("fft") lädt die Funktion.

Funktion: recttopolar (a, b)

Transformiert komplexe Zahlen der Form a + b %i in die Polarform r %e^(%i t). a ist der Realteil und b der Imaginärteil der komplexen Zahl. Die Argumente a und b sind eindimensionale Arrays derselben Größe. Die Größe der Arrays muss eine Potenz von 2 sein.

Die Werte der originalen Arrays werden durch den Betrag r = sqrt(a^2 + b^2 und die Phase t = atan2(b, a) ersetzt. Die Phase ist ein Winkel in dem Bereich -%pi bis %pi.

recttoploar ist die inverse Funktion zu polartorect. Das Kommando load("fft") lädt die Funktion.

Funktion: inverse_fft (y)

Berechnet die inverse schnelle Fourier-Transformation. Das Argument y ist eine Liste oder ein Array mit den Daten, die zu transformieren sind. Die Anzahl der Daten muss eine Potenz von 2 sein. Die Elemente müssen Zahlen (ganze, rationale, Gleitkommazahlen oder große Gleitkommazahlen) oder numerische Konstanten sein. Weiterhin können die Elemente komplexe Zahlen a + b*%i sein, wobei der Realteil und der Imaginärteil wiederum Zahlen oder numerische Konstanten sein müssen.

inverse_fft gibt ein neues Objekt vom selben Typ wie y zurück. Die Ergebnisse sind immer Gleitkommazahlen oder komplexe Zahlen a + %i*b, wobei a und b Gleitkommazahlen sind.

Die inverse diskrete Fourier-Transformation ist wie folgt definiert. Wenn x das Ergebnis der inversen Fourier-Transformation ist, dann gilt für j von 0 bis n-1

x[j] = sum(y[k] exp(2 %i %pi j k / n), k, 0, n - 1)

Mit dem Kommando load("fft") wird die Funktion geladen. Siehe auch fft für die schnelle Fourier-Transformation.

Beispiele:

Reelle Daten.

(%i1) load ("fft") $
(%i2) fpprintprec : 4 $
(%i3) L : [1, 2, 3, 4, -1, -2, -3, -4] $
(%i4) L1 : inverse_fft (L);
(%o4) [0.0, 14.49 %i - .8284, 0.0, 2.485 %i + 4.828, 0.0, 
                       4.828 - 2.485 %i, 0.0, - 14.49 %i - .8284]
(%i5) L2 : fft (L1);
(%o5) [1.0, 2.0 - 2.168L-19 %i, 3.0 - 7.525L-20 %i, 
4.0 - 4.256L-19 %i, - 1.0, 2.168L-19 %i - 2.0, 
7.525L-20 %i - 3.0, 4.256L-19 %i - 4.0]
(%i6) lmax (abs (L2 - L));
(%o6)                       3.545L-16

Komplexe Daten.

(%i1) load ("fft") $
(%i2) fpprintprec : 4 $                 
(%i3) L : [1, 1 + %i, 1 - %i, -1, -1, 1 - %i, 1 + %i, 1] $
(%i4) L1 : inverse_fft (L);
(%o4) [4.0, 2.711L-19 %i + 4.0, 2.0 %i - 2.0, 
- 2.828 %i - 2.828, 0.0, 5.421L-20 %i + 4.0, - 2.0 %i - 2.0, 
2.828 %i + 2.828]
(%i5) L2 : fft (L1);
(%o5) [4.066E-20 %i + 1.0, 1.0 %i + 1.0, 1.0 - 1.0 %i, 
1.55L-19 %i - 1.0, - 4.066E-20 %i - 1.0, 1.0 - 1.0 %i, 
1.0 %i + 1.0, 1.0 - 7.368L-20 %i]
(%i6) lmax (abs (L2 - L));                    
(%o6)                       6.841L-17
Funktion: fft (x)

Berechnet die schnelle Fourier-Transformation. Das Argument x ist eine Liste oder ein Array mit den Daten, die zu transformieren sind. Die Anzahl der Elemente muss eine Potenz von 2 sein. Die Elemente müssen Zahlen (ganze, rationale, Gleitkommazahlen oder große Gleitkommazahlen) oder numerische Konstanten sein. Weiterhin können die Elemente komplexe Zahlen a + b*%i sein, wobei der Realteil und der Imaginärteil wiederum Zahlen oder numerische Konstanten sein müssen.

inverse_fft gibt ein neues Objekt vom selben Typ wie x zurück. Die Ergebnisse sind immer Gleitkommazahlen oder komplexe Zahlen a + %i*b, wobei a und b Gleitkommazahlen sind.

Die diskrete Fourier-Transformation ist wie folgt definiert. Wenn y das Ergebnis der Fourier-Transformation ist, dann gilt für k von 0 bis n-1

y[k] = (1/n) sum(x[j] exp(-2 %i %pi j k / n), j, 0, n - 1)

Sind die Daten x reelle Zahlen, dann werden die reellen Koeffizienten a und b so berechnet, dass gilt

x[j] = sum (a[k] * cos (2*%pi*j*k / n) 
          + b[k] * sin (2*%pi*j*k / n), k, 0, n/2)

wobei

a[0] = realpart (y[0])
b[0] = 0

und für k von 1 bis n/2-1

a[k] = realpart (y[k] + y[n - k])
b[k] = imagpart (y[n - k] - y[k])

sowie

a[n/2] = realpart (y[n/2])
b[n/2] = 0

Das Kommando load("fft") lädt die Funktion. Siehe auch inverse_fft für die inverse schnelle Fourier-Transformation.

Beispiele:

Reelle Daten.

(%i1) load ("fft") $
(%i2) fpprintprec : 4 $
(%i3) L : [1, 2, 3, 4, -1, -2, -3, -4] $
(%i4) L1 : fft (L);
(%o4) [0.0, - 1.811 %i - .1036, 0.0, .6036 - .3107 %i, 0.0, 
                         .3107 %i + .6036, 0.0, 1.811 %i - .1036]
(%i5) L2 : inverse_fft (L1);
(%o5) [1.0, 2.168L-19 %i + 2.0, 7.525L-20 %i + 3.0, 
4.256L-19 %i + 4.0, - 1.0, - 2.168L-19 %i - 2.0, 
- 7.525L-20 %i - 3.0, - 4.256L-19 %i - 4.0]
(%i6) lmax (abs (L2 - L));
(%o6)                       3.545L-16

Komplexe Daten.

(%i1) load ("fft") $
(%i2) fpprintprec : 4 $
(%i3) L : [1, 1 + %i, 1 - %i, -1, -1, 1 - %i, 1 + %i, 1] $
(%i4) L1 : fft (L);
(%o4) [0.5, .3536 %i + .3536, - 0.25 %i - 0.25, 
0.5 - 6.776L-21 %i, 0.0, - .3536 %i - .3536, 0.25 %i - 0.25, 
0.5 - 3.388L-20 %i]
(%i5) L2 : inverse_fft (L1);
(%o5) [1.0 - 4.066E-20 %i, 1.0 %i + 1.0, 1.0 - 1.0 %i, 
- 1.008L-19 %i - 1.0, 4.066E-20 %i - 1.0, 1.0 - 1.0 %i, 
1.0 %i + 1.0, 1.947L-20 %i + 1.0]
(%i6) lmax (abs (L2 - L));
(%o6)                       6.83L-17

Berechnung der Sinus- und Kosinus-Koeffizienten.

(%i1) load ("fft") $
(%i2) fpprintprec : 4 $
(%i3) L : [1, 2, 3, 4, 5, 6, 7, 8] $
(%i4) n : length (L) $
(%i5) x : make_array (any, n) $
(%i6) fillarray (x, L) $
(%i7) y : fft (x) $
(%i8) a : make_array (any, n/2 + 1) $
(%i9) b : make_array (any, n/2 + 1) $
(%i10) a[0] : realpart (y[0]) $
(%i11) b[0] : 0 $
(%i12) for k : 1 thru n/2 - 1 do
   (a[k] : realpart (y[k] + y[n - k]),
    b[k] : imagpart (y[n - k] - y[k]));
(%o12)                        done
(%i13) a[n/2] : y[n/2] $
(%i14) b[n/2] : 0 $
(%i15) listarray (a);
(%o15)          [4.5, - 1.0, - 1.0, - 1.0, - 0.5]
(%i16) listarray (b);
(%o16)           [0, - 2.414, - 1.0, - .4142, 0]
(%i17) f(j) := sum (a[k] * cos (2*%pi*j*k / n) + b[k] 
                         * sin (2*%pi*j*k / n), k, 0, n/2) $
(%i18) makelist (float (f (j)), j, 0, n - 1);
(%o18)      [1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0]
Funktion: horner (expr, x)
Funktion: horner (expr)

Formt ein Polynom expr in das Horner-Schema um. Mit x wird die Variable angegeben, für die das Horner-Schema zu bilden ist. Wird das Argument x nicht angegeben, wird die Hauptvariable des kanonischen Ausdrucks expr für die Bildung des Horner-Schemas genutzt.

Das Horner-Schema kann die Stabilität der numerischen Berechnung eines Ausdrucks verbessern.

Beispiel:

(%i1) expr: 1e-155*x^2 - 5.5*x + 5.2e155;
                           2
(%o1)            1.0E-155 x  - 5.5 x + 5.2E+155
(%i2) expr2: horner (%, x), keepfloat: true;
(%o2)            (1.0E-155 x - 5.5) x + 5.2E+155
(%i3) ev (expr, x=1e155);
Maxima encountered a Lisp error:

 floating point overflow

Automatically continuing.
To reenable the Lisp debugger set *debugger-hook* to nil.
(%i4) ev (expr2, x=1e155);
(%o4)                       7.0E+154
Funktion: find_root (expr, x, a, b)
Funktion: find_root (f, a, b)
Optionsvariable: find_root_error
Optionsvariable: find_root_abs
Optionsvariable: find_root_rel

Findet die Nullstellen eines Ausdrucks expr oder einer Funktion f in dem Intervall [a, b]. Der Ausdruck expr kann eine Gleichung sein. In diesem Fall sucht die Funktion find_root die Nullstellen für den Ausdruck lhs(expr) - rhs(expr).

Kann Maxima den Ausdruck expr oder die Funktion f in dem Intervall [a, b] für alle Werte auswerten und ist der Ausdruck expr oder die Funktion f in dem Intervall stetig, dann ist sicher, dass find_root die Nullstelle oder zumindest eine Nullstelle findet, wenn mehrere Nullstellen vorhanden sind.

find_root beginnt mit einer binären Suche der Nullstelle. Erscheint die Funktion als glatt genug, wendet Maxima einen Algorithmus mit einer linearen Interpolation für die Suche der Nullstelle an.

Die Genauigkeit der Nullstellensuche wird von den Optionsvariablen find_root_abs und find_root_rel kontrolliert. find_root endet, wenn die Auswertung der Funktion ein Ergebnis hat, das kleiner als find_root_abs ist oder wenn aufeinander folgende Auswertungen Ergebnisse x_0 und x_1 haben, die sich voneinander weniger als find_root_rel * max(abs(x_0), abs(x_1)) unterscheiden. Der Standardwert der Optionsvariablen find_root_abs und find_root_rel ist Null.

find_root erwartet, dass die Funktion an den Endpunkten des Intervalls für die Nullstellensuche ein unterschiedliches Vorzeichen hat. Hat die Funktion an den Endpunkten des Intervalls dasselbe Vorzeichen, wird das Verhalten der Funktion find_root von der Optionsvariablen find_root_error kontrolliert. Hat find_root_error den Wert true, wird eine Fehlermeldung ausgegeben. Ansonsten wird von find_root der Wert von find_root_error als Ergebnis zurückgegeben. Der Standardwert von find_root_error ist true.

Kann die Funktion f bei der Nullstellensuche nicht zu einer Zahl ausgewertet werden, gibt find_root ein teilweise ausgewertetes Ergebnis zurück.

Die Reihenfolge der Grenzen des Intervalls a und b wird ignoriert. find_root sucht die Nullstellen immer in dem Intervall [min(a, b), max(a, b)].

Beispiele:

(%i1) f(x) := sin(x) - x/2;
                                        x
(%o1)                  f(x) := sin(x) - -
                                        2
(%i2) find_root (sin(x) - x/2, x, 0.1, %pi);
(%o2)                   1.895494267033981
(%i3) find_root (sin(x) = x/2, x, 0.1, %pi);
(%o3)                   1.895494267033981
(%i4) find_root (f(x), x, 0.1, %pi);
(%o4)                   1.895494267033981
(%i5) find_root (f, 0.1, %pi);
(%o5)                   1.895494267033981
(%i6) find_root (exp(x) = y, x, 0, 100);
                            x
(%o6)           find_root(%e  = y, x, 0.0, 100.0)
(%i7) find_root (exp(x) = y, x, 0, 100), y = 10;
(%o7)                   2.302585092994046
(%i8) log (10.0);
(%o8)                   2.302585092994046

Nächste: , Vorige: , Nach oben: Fourier-Transformationen   [Inhalt][Index]