1 − 2 53 liegt. Desweiteren

Warum selber rehnen es gibt doh Computer?

Teilnehmer:

Christian Draeger Andreas-Obershule

Dennis Menge Heinrih-Hertz-Obershule

Patrik Mikut Immanuel-Kant-Obershule

Andre Rozek Herder-Obershule

Patrik Klitzke Andreas-Obershule

ChristopherShmidt Immanuel-Kant-Obershule

Robert Shlämann Immanuel-Kant-Obershule

Gruppenleiter:

FalkEbert Tehnishe Universität Berlin,

Mitglied imDFG-Forshungszentrum Matheon

Mathematikfür Shlüsseltehnologien

Wenn man einmal die Shüler eines Mathematikkurses beobahtet, dann stellt

man shnell fest, dass sie ihre Tashenrehner öfter benutzen als ihren Kopf.

Wozu man sih auh die Mühe mahen, einen Wert wie

(2 ⁵⁰ + 1) − 2 ⁵⁰

^selbst

auszurehnen, wenn der Tashenrehner shon gribereit liegt?

Warum man seinem Tashenrehner nihtimmer trauenkann, haben wir in die-

semeinwöhigenKurserfahren.Dazumuss manersteinmalwissen,wieeinReh-

ner überhaupt rehnet, weshalb wir uns zuerst mitder Mashinenarithmetikbe-

fasst haben. Weiterhin haben wir die möglihen Rehenoperationen betrahtet,

dennfür einen Computer sind Addition und Multiplikation mitrelativgeringem

Aufwand möglih.Wieman aber dieDivision und dieExponentialfunktion oder

auhtrigonometrisheFunktionennuraufdiesebeidenOperationenzurükführen

kann, wollen wir in unserem Beriht zeigen. Um dies zu ermöglihen, benutzten

wir Approximationen mit dem Newton-Verfahren und mit Taylor-Polynomen.

Neben der Möglihkeit, Werte selbst zu errehnen, haben wir gezeigt, wie es

durh Interpolation von Tabellen möglih ist, alle Werte exakt und ezient zu

bestimmen.

2 Mashinenarithmetik

Esgibtvershiedene Zahlensysteme.Dasgängistefürunsistdas Dezimalsystem.

IndiesemSystem stehen uns dieZiern 0bis 9zur Verfügung. EinemComputer

bzw.einemprozessorgesteuertenGerätstehenabernurzweiZuständezur Verfü-

gung,nämlihanundaus.DieseZuständewerdenübliherweisedurhdieZiern

0(füraus) und 1(füran) dargestellt.Das Zahlensystem,das nuraus den Ziern

0und 1 bestehtwird alsDual-oder auhBinärsystem bezeihnet.

2.1 Darstellung von Dualzahlen

Da die Darstellung und Umrehnung von Ganzzahlen im Dualsystem allgemein

bekannt ist, möhten wir in unserem Beriht nur auf die Darstellung von Flieÿ-

kommazahleneingehen.Eine typishe normierteFlieÿkommazahlwirdfolgender-

maÿenbeshrieben:

± 0, 1m 2 m 3 ..m 53 · 2 ^{e 10 e 9 ..e 0 − 1023}

Hierbeihandeltessihumein64bit-System.DasersteBitwirdfürdasVorzeihen

derDualzahlverwendet:Eine0bedeutet,dassdieZahlpositivist,eine1dagegen,

dass sie negativ ist. Die Mantisse der Dualzahl beginnt mit 0,1. Dies ist eine

Festlegung, um die Eindeutigkeit der Zahlen zu garantieren. Danah folgen 52

Ziern der Mantisse, welhe im Bereih von

0, 5

^bis

1 − 2 ⁵³

^liegt. Desweiteren

Zahlwird dann 1023 subtrahiert.Der Exponentist danneine Ganzzahl, diesih

imBereih von 1024 bis-1023 bendet.

Desweiteren gibt esfolgende Festlegungen:

m ₂ ..m ₅₃ = 0 ∧ E = 0 → 0

⁽¹⁾

m ₂ ..m ₅₃ = 0 ∧ E = 2047 → ±∞

⁽²⁾

mindestens ein

m i 6 = 0 ∧ E = 2047 → NaN (NotaNumber)

⁽³⁾

E

^{beizeihnet hierbei die durh}

e 0

^bis

e 10

^{gebilete Zahl und NaN eine Zahl, die}

gar keine ist,wie z.B.

0 0

^.

ImweiterenVerlaufdiesesKapitels,wirddas8Bit-Systemverwendet,wodieMan-

tisse5StellenundderExponent2Stellenbesitzt.Damitlassensihbeispielsweise

positiveZahlenimBereih von

0, 10000 · 2 ^{− 11} = 0, 0625

^bis

0, 11111 · 2 ⁺¹¹ = 7, 75

darstellen.

2.2 Mashinengenauigkeit

JederComputermahtbeimRehnenFehler.Mathematishbezeihnet mandie-

sen Fehler wie folgt:

eps = b 2 · b ^{− µ} ,

wobei

µ

^{die Anzahl der Stellender Mantisse und b dieBasis des} Zahlensystems angibt.Speziell für das Dualsystem gilt:

eps = 2

2 · 2 ^{− 5} = 1 32

Im8 Bit System istdieMashinengenauigkeitalso etwasmehr alseine Dezimal-

stelle.Allgemein giltFolgendes:

1.

rd(x) = x(1 + ε 1 )

2.

x ⊗ y = (x × y)(1 + ε ₂ ), × ∈ { +, − , · , ÷}

3.

| ε 1/2 | ≤ eps

Dies heiÿt, dass es sowohl zu Fehlern bei der Darstellung von Dezimalzahlen in

dualen Zahlen kommt, als auh bei jeder möglihen Kombinationvon mehreren

Dualzahlen.Dieserelativen Fehlerkönnenmaximalso groÿsein, wie dieMashi-

nengenauigkeit.

Umwandlung von

π

^{von der dezimalen Zahl in die duale Zahl durh das 8Bit-}

Systemdemonstrieren. DurhdieMashiengenauigkeitvon

1

32

^{ergibt sih,dass}

π

zwishen

31

32 π ≈ 3, 04

^und

³³ ₃₂ π ≈ 3, 24

^liegt.

Zuerstberehnen wir den Exponenten. Eswird bald ersihtlih, dass

π

^zwishen

0, 5 · 2 ² = 2

^und

0, 5 · 2 ³ = 4

^{liegt . Demnahist der Exponent 2. Mankann dies}

auh durh den Logarithmusberehnen:

0.5 · 2 ^x = π log(0.5 · 2 ^x ) = log(π) log(0.5) + x · log(2) = log(π)

x = log(π) − log 0.5 log(2) x ≈ = 2.65

Hierbeiwird abgerundet, dadieeigentlihe Mantisse gröÿer als0,5 sein kann.

Nun muss nur noh die Mantisse berehnet werden. Dazu teilen wir die Zahl

durh

2 ^Exponent

^{. Wir erhalten rund 0,7854. Nun wird diese Zahl mit zwei multi-}

pliziert. Wenn Sie gröÿer als eins ist, ist die nähste Stelle der Mantisse 1

eins

und eins wird abgezogen, anderenfalls null. Dies wird nun solange gemaht, bis

alleMantissenstellen gefülltsind:

π : 2 ² ≈ 0, 7854 0, 7854 · 2 = 0, 5708 + 1 0, 5708 · 2 = 0, 1416 + 1 0, 1416 · 2 = 0, 2832 + 0 0, 2832 · 2 = 0, 5664 + 0 0, 5664 · 2 = 0, 1328 + 1 0, 1328 · 2 = 0, 2656 + 0 0, 11001 | 0

Dadasarry-Bitnullist,mussnihtgerundetundauhnihtnormalisiertwerden.

Wirerhalten nun

0, 11001 · 2 ¹⁰ = 0, 7854 · 2 ² = 3.125

^{.Wir erhalten also3.125 im}

8BitSystem für

π

^.

2.3 Addition mit Dualzahlen

Vorgehensweise Um zwei Dualzahlen miteinander addieren zu können, müs-

sensiedengleihenExponentenbesitzen.Derkleinerewirdandengröÿerenange-

1

InunseremFalldieerste

(also das Komma nah links vershoben) wird. Dieser Vorgang wird 'Shiften'

genanntundistauhindieandereRihtungmöglih.Addiertmannun dieMan-

tissen, muss die Dualzahl eventuell normalisiert werden, denn durh Überträge

kann diese von der üblihen Darstellungsweise abweihen (z.B.

1, m ₂ m ₃

^..).Dies

geshiehtebenfallsdurh eineadditiveVeränderung desExponenten beientspre-

hender Multiplikation bzw. Division der Mantisse mit2. Anshlieÿend wird die

Mantisse gerundet, was mitHilfe des sogenannten 'arry-Bits' geshieht.

Dasarry-Bit wirdnurkurzzeitiggespeihert und gibt die6.Stelle der Mantisse

an. Beträgt das arry-Bit 0, so wird abgerundet, anderenfalls aufgerundet. Die

Aufrundungkann dazu führen, dass dieMantisse gleih 1 wird und einweiteres

malnormalisiertwerden muss.

Beispielrehnung Nehmen wir nun ein einfahes Beispiel: Gegeben seien die

zwei Zahlen3und 1,4375.Wenn wirdiese inBinärzahlenumwandeln,ergibtsih

0, 11000 · 2 ¹⁰

^und

0, 10111 · 2 ⁰¹

^{.DiekleinereZahlwirdnunnahrehtsvershoben:}

0, 010111 · 2 ¹⁰

^{. Nunaddieren wir diese:}

0, 11000 2 ¹⁰ +0, 010111 2 ¹⁰ 1, 00011 | 1 2 ¹⁰

DanundieMantissegröÿeralseinsist,wirdsienormaliesiert.Wirvershiebendie

Mantisseum eins nahrehts 2

und erhöhenden Exponenten um eins:

0, 1000111 · 2 ¹¹

^{.Da das arry-Bit3 1ist, wird}aufgerundet. Daraus ergibt sih:

0, 10010 · 2 ¹¹

^.

DadieMantisse nihtgleiheins ist,mussniht noheinmalnormalisiertwerden.

Wir erhalten

0, 10010 · 2 ¹¹ = 4, 5

^{. Der rihtige Wert beträgt 4,4375. Der Fehler}

entsteht aufgrund der Mashinengenauigkeit.

2.4 Multiplikation von Dualzahlen

Vorgehensweise Die Dualzahlen werden multipliziert,indemdie Exponenten

miteinanderaddiertunddieMantissenmiteinandermultipliziertwerden.Wiebei

der Addition wird die Dualzahl nun miteiner geeigneten Shiftoperation norma-

lisiert,anshlieÿendgerundet und ggf.ein weiteresmal normalisiert.

2

Division durhzwei

3

6.Stelle

dieAnshalihkeitund dem besseren Verständnis sehr hilfreih.

Aus Gründen der Einfahheit nehmen wir wieder die gleihen Zahlen wie auh

beimAdditionsbeispiel.Wir multiplizieren3 mit1,4375 im8Bit-System:

0, 10111 · 0, 11000 = 0, 10111 · 0, 1 + 0, 10111 · 0, 01

= 0, 10001 | 01

Jetzt muss gerundet werden: Dabei wird die sehste und siebte Stelle der Man-

tisse abgeshnitten. Es wird niht aufgerundet, da das arry-Bit null ist. An-

shlieÿend müssen die Exponenten noh addiert werden. Es ergibt sih dann

0, 10001 · 2 ¹⁰⁺⁰¹ = 0, 10001 · 2 ¹¹ = 0.53125 · 2 ³ = 4.25

^{.Dieserstimmtgerundet mit}

der ersten Stelledes tatsählihen Wertes von 4,3125 überein.

3 Approximation von Funktionswerten

WirwerdenimFolgenden einigeTehniken zumnäherungsweisen Auswertenvon

analytishen Funktionen entwikeln. Die Standard-Methode dazu ist dieTaylor-

entwiklung. Andere spezielle Funktionswerte - insbesondere die Werte von Um-

kehrfunktionen - lassensihmitHilfe vonNullstellensuhbestimmen.Zudiesem

Zwek wird das Newton-Verfahren vorgestellt.

3.1 Taylor-Entwiklung

Wir betrahten eine genügend oft dierenzierbare Funktion

f

^{in der Nähe des}

Punktes

x 0

^{. Es gilt dann, wenn man genügend nahe an}

x 0

^{bleibt, dass}

f (x) ≈ f(x 0 )

^{. Dabei setzen wir}

p(x) = f (x 0 )

^{. Dies ist eine konstante Funktion. Eine}

bessere Näherung an

f

^{in der Nähe von}

x 0

^{ist zu erwarten, wenn man statt}

einer konstanten eine lineare Funktion

p 1

^{zulässt. Dabei fordern wir zusätzlih}

zu

p ₁ (x ₀ ) = f (x ₀ )

^{, dass}

p ^′ ₁ (x ₀ ) = f ^′ (x ₀ )

^{. Diese Gerade ist die Tangente an die}

Funktion

f

^imPunkt

x 0

^{.Eine Gleihung für}

p 1

^{ergibtsihmitder} Grenzwertbe- trahtung

x lim → x ⁰

f (x) − p 0 (x) x − x 0

.

⁽⁴⁾

Zählerund NennerdiesesAusdruks konvergieren gegen0,da

p 0 (x) = f (x 0 )

^.Mit

dem Satz von L'Hospital erhältman

x lim → x 0

f (x) − p 0 (x)

x − x ₀ = lim

x → x 0

f ^′ (x) − 0

1 − 0 = f ^′ (x 0 ).

Ausmultiplizierenliefert

x lim → x ⁰ f(x) = lim

x → x ⁰ f (x 0 ) + f ^′ (x 0 )(x − x 0 )

^def

= lim

x → x ⁰ p 1 (x).

Diese Gerade stimmt in

x 0

^mit

f

^{sowohl im Wert als auh in der ersten Ablei-}

tung überein. Wir versuhen jetzt eine quadratishe Parabel zu nden, die mit

f

^{inWert,erster und zweiter Ableitung}übereinstimmt. Dazubetrahten wir die Dierenz

f(x) − p 1 (x)

^{.An der Stelle}

x 0

^{isdiese Dierenz nullund auhdieerste}

Ableitung.

f (x) − p ₁ (x)

^{hat alsoeinenkritishen Punkt an}

x ₀

^{und verhältsih in}

einerkleinenUmgebung wie

α(x − x 0 ) ²

^.UmdieZahl

α

^{zu bestimmen,betrahten}

wir

x lim → x 0

f (x) − p 1 (x) (x − x 0 ) ² .

Diemehrfahe Anwendung des Satzes von L'Hospital erhältman

x lim → x ⁰

f(x) − p 1 (x)

(x − x 0 ) ² = lim

x → x ⁰

f(x) − (f(x 0 ) + f ^′ (x 0 )(x − x 0 ))

(x − x 0 ) ²

⁽⁵⁾

= lim

x → x 0

f ^′ (x) − (f ^′ (x 0 )) 2(x − x 0 )

= lim

x → x 0

f ^′′ (x)

2 = f ^′′ (x 0 ) 2 .

Ausmultiplizierenliefert

x lim → x 0

f (x) = lim

x → x 0

f (x 0 ) + f ^′ (x 0 )(x − x 0 ) + 1

2 f ^′′ (x 0 )(x − x 0 ) ²

^def

= lim

x → x 0

p 2 (x).

Mitder analogen Vorgehensweise wie in (4) und (5) erhältman für beliebige

n p _n (x) = f(x ₀ )

0! + f ^′ (x ₀ )

1! (x − x ₀ ) + f ^′′ (x ₀ )

2! (x − x ₀ ) ² + . . . + f ⁽ⁿ⁾ (x ₀ )

n! (x − x ₀ ) ⁿ .

⁽⁶⁾

DieseNäherung an

f

^{inder Nähe von}

x 0

^{nennt sih T}aylor-Polynom.

DieGenauigkeitder Approximation von

f

^durh

p n

^{zeigenwir imFolgenden.Es}

gilt

(f(x) − p n (x)) ⁽ⁿ⁺¹⁾ = f ⁽ⁿ⁺¹⁾ (x) − p ⁽ⁿ⁺¹⁾ _n = 0

weil

p n

^einPolynom

n

^{-ten Gradesist unddie}

(n + 1)

^{teAbleitunggleihnull ist.}

Weiterhin gilt

(f (x) − p n (x)) ⁽ⁿ⁾ = Z x

x 0

(f (x) − p n (x)) ⁽ⁿ⁺¹⁾ = (f ⁽ⁿ⁾ (x) − f ⁽ⁿ⁾ (x 0 )) − (p ⁽ⁿ⁾ _n (x) − p ⁽ⁿ⁾ _n (x 0 )).

Mitdem Mittelwertsatzder Dierentialrehnung ergibt sih

(f (x) − p n (x)) ⁽ⁿ⁾ = f ⁽ⁿ⁺¹⁾ (ξ)(x − x 0 ) − p ⁽ⁿ⁺¹⁾ _n (ξ)

| {z }

= 0

⁽⁷⁾

fürein

ξ ∈ [x, x 0 ]

^{. Die}

n

^maligeIntegration von(7) liefertden Ausdruk

f(x) − p n (x) = f ⁽ⁿ⁺¹⁾ (ξ)

(n + 1)! (x − x 0 ) ⁿ⁺¹ .

Wirnehmen

f(x) = e ^x

^und

x ₀ = 0

^{. Dannergibt sih mit(6)}

p n (x) = 1 + x + x

2 + x

6 + . . . + x ⁿ

n! .

⁽⁸⁾

Undes gilt

e ^x − p n (x) = e ^ξ x ⁿ⁺¹

(n + 1)! ≤ max(1, e ^ξ ) x ⁿ⁺¹ (n + 1)! .

Um

e ^x

^mit

x = M · 2 ^E

^{zu bestimmen, kann man wiefolgt vorgehen. Es gilt}

| M | ∈ [ ¹ ₂ , 1)

^und

E ∈ Z

.

e ^x = (e ^M ) ^{2 E} .

⁽⁹⁾

Esgenügtalso,dieExponentialfunktionfürWertein

[ − 1, 1]

^{auswertenzukönnen.}

Dies kann mit Hilfe der Taylorpolynome

p _n

^{geshehen. Dabei ist nah (7) der}

Fehler maximal

e

(n+1)!

^{. Für 8 Stellen} Genauigkeit benötigt man shon

p 12

^{. Die}

Nutzungder Potenzgesetze spart aberRehenzeit.

e ^M = (e ^M/256 ) ²⁵⁶ = ((((((((e ^M/256 ) ² ) ² ) ² ) ² ) ² ) ² ) ² ) ² .

⁽¹⁰⁾

ˆ Bestimme

M 2 ⁸

^.

ˆ Bestimme

p 4 (M )

^mit

p 4

^{aus (8).}

ˆ Quadriere

p ₄ (M ) | E | + 8

^mal.

ˆ Wenn

E

^{negativ ist,bestimmedas Reziprokedes letzten}Ergebnisses.

ˆ Passe das Ergebnisandie Rehnerarithmetikan.

Dabei sind ein8-StellenShift sowie 4 Additionen und 4Multiplikationenfür die

Auswertung von

p 4

^{nötig. Hinzu kommen weitere}

| E | + 8

Multiplikationen und eine potentielle Division. Da selbst in einer 64bit Arithmetik

e ^{2 10} = Inf

^{ist, ist}

auh dieAnzahl an Multiplikationenauf 22beshränkt.

Fürhöhere Genauigkeiten kann das benutzte

p n

^{nohweiter erhöht werden.}

3.2 Newton-Verfahren

Das Finden vonNullstellen von Funktionenist niht immereinfah. Ausdiesem

Grundgibt es Nährungsverfahren, wie z.B. das Newton-Verfahren.

Das Newton-Verfahren funktioniert am besten, wenn als Startpunkt ein Punkt

gewählt wird, der möglihst nah ander gesuhten Nullstelle liegt.Wirversuhen

in diesem Verfahren, die Nullstelle

x ^∗

^{durh mehrere} Iterationen anzunähern.

a

³

M 0, 75

E 2

M/256 1, 4142136 p 4 (M/256) 1, 00293398 (p 4 (M/256)) ² 56 2, 11700002 ((p 4 (M/256)) ² 56) ⁴ 20, 0855369 M 0, 627673039

E 5

e ³ 20, 0855369

Tabelle1:Shritte beider Bestimmung von

e ³

Voraussetzungen dafür, dass das Newton-Verfahren überhaupt funktioniert,sind

dieFolgenden:

ˆ eine stetig dierenzierbareFunktion, diemindestens eine Nullstelle besitzt

ˆ dieerste Ableitungsollte niht 0 werden ander Nullstelle

Esgibt vershiedene Herleitungenfür dieses Verfahren. ImFolgenden stellen wir

einevondiesenvor.DiesebasiertaufderzuvorvorgestelltenTaylor-Entwiklung.

Sei

f(x)

^{eine stetig} dierenzierbare Funktion. Weiterhin sei

f(x ^∗ ) = 0

^{, al die}

Funktion

f

^{besitzt eine Nullstellebei}

x ^∗

^{. Dannsei}

x 0 ≈ x ^∗

^{.Jetzt betrahten wir}

das

p 1

^{-Polynom der T}aylorentwiklung, alsodas Polynom biszum Glied mitder 1.Ableitung.

f(x) ≈ f(x 0 ) + f ^′ (x 0 ) · x − x 0 = p 1 (x)

Jetzt bestimmen wir die Nullstelle

x 1

^{dieses Polynoms.}

f(x ₀ ) + f ^′ (x ₀ ) · (x ₁ − x ₀ ) = 0

Nunstellt man dieGleihung einfah nur noh nah

x 1

^um.

f ^′ (x ₀ ) · (x ₁ − x ₀ ) = − f (x ₀ ) x ₁ − x ₀ = − f(x 0 ) f(x 0 ) x 1 = x 0 − f (x 0 )

f ^′ (x 0 )

herungenfür

x ^∗

^{zu erhalten undesergibtsihdannfolgende}Iterationsvorshrift:

x n+1 = x n − f (x n )

f ^′ (x n )

⁽¹¹⁾

3.3 Konvergenz des Newton-Verfahren

Wirbetrahten dieIterationsfunktion des Newton-Verfahrens

ϕ(x) = x − f (x) f ^′ (x) .

Esgilt fürdie Nullstelle

x ^∗

^dass

ϕ(x ^∗ ) = x ^∗

^{. und weiterhin}

ϕ ^′ (x) = 1 − (f ^′ (x)) ² − f(x) · f ^′′ (x)

(f ^′ (x)) ²

= 1 − (f ^′ (x)) ²

(f ^′ (x)) ² + f (x) · f ^′′ (x) (f ^′ (x)) ²

= f(x) · f ^′′ (x)

(f ^′ (x)) ² .

⁽¹²⁾

Damitgilt

x ^∗ − x _n+1 = x ^∗ − ϕ(x _n )

= ϕ(x ^∗ ) − ϕ(x n )

= ϕ(x ^∗ ) − (ϕ(x ^∗ ) + ϕ ^′ (x ^∗ )(x n − x ^∗ ) + ϕ(ξ)

2 (x n − x ^∗ ) ² .

⁽¹³⁾

Dabei ist

ϕ(x n )

^{bis zum Polynom}

p 1

^{der T}aylorentwiklung um

x ^∗

^entwiklelt

und das Fehlerglied addiertworden. Mit (12) folgt

ϕ ^′ (x ^∗ ) = 0

^da

f (x ^∗ ) = 0

^und

f ^′ (x ^∗ ) 6 = 0

^{.Für (12) folgt damit}

x n+1 − x ^∗ = ϕ ^′′ (ξ)

2 (x n − x ^∗ ) ²

für ein

ξ ∈ [x ^∗ , x n ]

^{. Da der Abstand zur Nullstelle im}

(n + 1)

^{ten Shritt vom}

Quadratdes Abstandsdes

n

^{tenShrittesabhängt, spriht manauhvonquadra-}

tisherKonvergenz.Wenn

(x n − x ^∗ )

^{genügendkleinist,dannist}

(x n+1 − x ^∗ )

^noh

sehr viel kleiner. Anshaulih kann man sagen, dass sih die Anzahl korrekter

Stellender Iterierten

x n

^{in jedem Shritt etwa} verdoppelt.

Die Division mit Hilfe des Newton-Verfahrens beruht auf dem Prinzip, dass wir

unseineFunktionssharsuhendiealsNullstellendieReziprokenwertederZahlen

a

^{hat, durh die wir dividieren möhten. Eine spezielle Wahl dieser Funktionen}

sorgt dafür, dass das Newton-Verfahren ohneDivisionen anwendbarist.

Zur Bestimmung des Reziproken einer Zahl

a

^{kann man wie folgt vorgehen. Zur}

Vereinfahung gehenwir von

a > 0

^{aus. Das Vorzeihen von}

a

^und

a ^{− 1}

^{ist iden-}

tish.Wir mahen diefolgendenVorbetrahtungen

ˆ

a = M · 2 ^E

^,

ˆ

1

a = _M ¹ · 2 ^{− E}

^,

ˆ

f(x) = ¹ _x − M

^hat

_M ¹

^alsNullstelle,

ˆ

x n+1 = x n −

1 x n − M

− ¹

x 2 n

= x n (2 − M · x n )

^.

Manmuss alsonur das Reziproke vonZahlen

M ∈ [0, 5; 1)

^{bestimmen können.}

a = M · 2 ^E x _n+1 = x _n −

1 x n − M

− _x ^{1 2} _n = x _n (2 − M · x _n )

DieseGleihungdientunszur Reziprokenberehnung.Dadas Verfahrenshneller

konvergiert,wennmanStartwertewählt,dieinder NähedergesuhtenNullstelle

liegen,benötigen wireinVerfahren,das unsgeeigneteStartwerte

(x 0 )

^liefert,also

solhe, die bereits nahe an

1

M

^{liegen. Eine einfahe Möglihkeit ist eine Aproxi-}

mationvon

1/M

^{durh eine lineareFunktion der Form}

g(x) = mx + n

^.

DieFunktionsgleihung kann man durh folgendesGleihunggssystem gewinnen:

f (x) = 1 x f( 1

2 ) − g( 1

2 ) = δ

f (1) − g(1) = δ

f (ξ) − g(ξ) = − δ

(f (ξ) − g(ξ)) ^′ = 0

0 0.5 1 1.5 2 2.5 0

0.5 1 1.5 2 2.5

x

y

1/x mit Minimax−Gerade

1/x g(x)

Abbildung1: Minimax-Gerade

Dabeifordernwir,dassdermaximaleAbstandzwishen

f

^und

g

^{anden Rändern}

und an einer Zwishenstelle angenommen wird und insgesamt minimal ist. Die

Lösungdieses Gleihungsystems istdie sogenannte Minimax-Gerade

g(x) = − 2 · x + 3 2 + √

2.

Mitdiesem Startwert reihen für 8Stellen Genauigkeit etwa3 Newton-Shritte.

ˆ Bestimme

x 0 = g(M ) = − 2 · M + ( ³ ₂ + √ 2)

^.

ˆ Bestimme

x 3 ≈ _M ¹

^.

ˆ Rehenbedarf: 7Multiplikationen,4 Additionen

Da

1/M ∈ (1; 2]

^{liegtmüssenwirnunnohdieZahlanden}Mantissenbereihund den Exponenten anpassen.

DiesesVerfahrenzurDivision(nahMultiplikationmitdementsprehenden Zäh-

ler) ist tatsählih shneller als die Vorgehensweise wie beim shriftlihen Divi-

dieren.

1

a ≈ x 3

2 · 2 ^{− E +1}

a

³

M 0, 75

E 2

x 0 1, 4142136 x 1 1, 3284271 x 2 1, 3333153 x 3 1, 3333333 M 0, 6666667

E − 1

1/a 0, 3333333

Tabelle2: Shritte beider Bestimmungvon

1/3

Anwendungen

WeitereAnwendungsmöglihkeitensindzum Beispieldas Logarithmieren.Esgilt

ˆ

a = M · 2 ^E

^,

ˆ

ln(a) = ln(M ) + ln(2) · E

^,

ˆ

e ^x − M

^hat

ln(M )

^alsNullstelle.

ˆ

x n+1 = x n − ^{e xn} e ^{xn − M} = x n − 1 + M · e ^{− x n}

^(Newton)

Es muss also nur

ln(M )

^für

M ∈ [ ¹ ₂ , 1)

^{durh einige} Newtoniterationen aus- gewertet werden. Dafür sind nur Additionen und Multiplikationen notwendig.

Zusätzlihmuss

ln(2)

^{alsKonstante zur Verfügungstehen. Die} Minimax-Gerade fürden Startwert ist

g(x) = 1, 38629436x − 1, 35646431

^.

4 Der CORDIC-Algorithmus

DerCORDIC 4

AlgorithmusisteineektivesVerfahren,um denCosinusundden

Sinus und damit den Tangens eines Winkels nährungsweise zu bestimmen. Der

Algorithmusbasiert auf folgenderIteration:

4

Die Abkürzung steht für Coordinate RotationDigital Computer und wurdevon Jak E.

VolderEndeder50erJahreentwikelt.

a

³

M 0, 75

E 2

x 0 − 0, 316743539 x 1 − 0, 287255667 x 2 − 0, 287681982 x 3 − 0, 287682073 ln(2) · 2 1, 38629436 x 3 + ln(2) · 2 1, 09861229 M 0, 549306144

E 1

ln(3) 1, 09861229

Tabelle3: Shritte beider Bestimmung von

ln3



 

 

 

 

x _n+1 = x _n − d _n · y _n · 2 ^{− n} y _n+1 = y _n + d _n · x _n · 2 ^{− n} z n+1 = z n − d n · arctan 2 ^{− n} d _n = sign(z _n )

(14)

Alle möglihen Werte von

arctan 2 ^{− n}

^{für die} vershiedenen

n

^{werden vorgespei-}

hert, wobei übliherweise

n

^{etwas gröÿer als die} Mantissenlänge gewählt wird.

Wenn

| z 0 |

^{kleiner oder gleih}

X ∞

k=0

arctan 2 ^{− k} = 1.743286047...

⁽¹⁵⁾

ist,dann ergibt sih

lim n →∞



 x n

y n

z n



 = K ×





x 0 · cos z 0 − y 0 · sin z 0

x 0 · sin z 0 + y 0 · cos z 0

0



 ,

wobei der Skalierungsfaktor K gleih

Π ^∞ _n=0 √

1 + 2 ^{− 2n} = 1.64676025...

^{ist. Um}

nun den Cosinus und Sinuseines Winkelszu bestimmen, benötigt man noh die

Startwerte für

x, y

^und

z

^{.Die Zahl für dieCosinus bzw. Sinusbestimmtwerden}

sollen,ist

θ

^.Dabeikann

θ

^{maximalsogroÿwiederbei(15)berehneteWertsein.}

x ₀ = 1

K = 0, 60725...

y 0 = 0 z ₀ = θ

DasGrundprinzipist,dassmanversuhtdurhDrehungeinesEinheitsvektorsum

vershiedene immerkleiner werdende Winkel

± arctan(2 ^{− n} )

einzukreisen.Ersetzt mandie Drehung um

± arctan(2 ^{− n} )

^durheine Drehstrekung mitdem zusätzli- henFaktor

√ 1 + 2 ^{− 2n}

^{,dannsinddiese} Drehstrekungen durhden Algorithmus

(14) realisiert. Zum Ausgleih der Strekung wird mit dem auf

1/K

^verkürzten

Vektor begonnen. Der nale Vektor ist der um

θ

^gedrehte Einheitsvektor, von dessenKomponenten

cos θ

^und

sin θ

^{abgelesen werden können.}

Die Vorteile bei diesem Verfahren sind zum einen, dass es sowohl Cosinus als

auhSinusgleihzeitigberehnetunddass manzumBerehnennuraddierenund

Mantissen shiften können muss.

5 Oder doh lieber Tabellen?

5.1 Grundlagen

Neben derMöglihkeitdieWerteeinerFunktion auszurehenen,istesauhmög-

lih, Tabellen zur Wertebestimmung heranzuziehen. Dabei werden bestimmte

Stellen die sogenannten Stützstellen und die dazugehörigen Werte gespei-

hert.DiesesVerfahrenwurdebevoresTashenrehner gab zur Bestimmungvon

beispielsweise Sinuswerten herangezogen. Damals ergab sih meist eine Genau-

igkeit von a. 4 Stellen,weil dietabellarisierten Werte niht genauer angegeben

waren. Um die Werte zwishen zwei Stützstellen zu bestimmen, wurde einfah

eine Gerade durh diese beiden Punkte gelegt und deren Werte als Näherungen

angenommen. Dies ist jedoh bei einer geforderten Mashinengenauigkeit von

mindestens 16 Stellen zu ungenau. Und besser als eine solhe Gerade wäre ein

Polynom

p(x) (n − 1)

^{-ten Grades, wobei}

n

^{dieAnzahl der} Stützstellen ist.

p(x) = a 0 + a 1 x + a 2 x ² + . . . + a n x ⁿ

Dieses Polynom muss an allen Stützstellen die in der Tabelle abgespeiherten

Werte annehmen,

p(x i ) = y i

^{. Abernoh istunklar obdieses Polynom überhaupt}

existiert.

Um die Existenz dieses Polynoms zu beweisen, denieren wir zuerst ein soge-

nanntes Lagrange-Polynom

(n − 1)

^{-Grades,für das gilt:}

L _i (x _j ) =

( 1, i = j

0, i 6 = j

⁽¹⁶⁾

Durh Multiplikation der Linearfaktorenergibt sih:

l i = (x − x 1 )(x − x 2 ) · · · (x − x i − 1 )(x − x i+1 ) · · · (x − x n ) l i (x i ) = (x i − x 1 )(x i − x 2 ) · · · (x i − x i − 1 )(x i − x i+1 ) · · · (x i − x n )

Dadurh ist in

l i

^{die Bedingung}

l i (x j ) = 0

^für

i 6 = j

^{erfüllt. Durh Teilung von}

l i

durh

l i (x i )

^{ist auh die zweite Bedingung}

l i (x)/l i (x i ) = 1

^für

i = j

^{erfüllt und}

L i (x)

^{ist nun deniert mit}

L i (x) = l i

l i (x i ) .

Daraus ergibt sihfür

p

p = a 1 L 1 + a 2 L 2 + . . . + a n L n .

Undfür dieStützstellen

x i

p(x i ) = a 1 L 1 (x i ) + a 2 L 2 (x i ) + . . . + a n L n (x i ).

Dabei ist jedoh wegen Bedingung (16) jeder Summand auÿer

a i L i (x i )

^gleih

Null, denn

L i

^und

a i

^{beziehen sih nur bei der} Stützstelle

x i

^{auf das gleihe i}

und damit ist der Faktor

L i

^{in jedem anderem Summanten Null. Folglih ist}

p(x n ) = a n = y n

^{. Durh} Koezientenvergleih an allenStützstellen ergibt sih

p(x) = y 1 L 1 (x) + y 2 L 2 (x) + . . . + y n L n (x).

DiesesPolynomexistiertunderfülltallevorheran

p(x)

^gestelltenAnforderungen.

Auh wenn die Existenz dieses Polynoms nun bewiesen ist, ist noh niht klar,

dass dieses eindeutig ist,dies wie folgtzeigen.

Nehmen wir einmal an, es gibt neben dem Polynom

p(x)

^{das Polynom}

q(x)

^,

welhes auh an allen Stützstellen

x i

^{durh die} zugehörigen Funktionswerte

y i

läuft.

q(x)

^istwie

p(x)

^einPolynom

(n − 1)

^{-ten Grades.}

Wären

p(x)

^und

q(x)

untershiedlih, dann lieÿe sih die Dierenz

r

^durh

r = p − q

beshreiben. Da dies anallen Stellen von

p

^und

q

^{gilt,folgt darausauh}

r(x _i ) = p(x _i ) − q(x _i ).

⁽¹⁷⁾

r(x)

^{isteinPolynom}

(n − 1)

^{-ten Grades,daesdieDierenzvonzweiPolynomen}

(n − 1)

^{-ten Grades ist. Aus (17) geht hervor, dass}

r(x i )

^{immer gleih Null ist,}

da

p

^und

q

^{anden Stellen}

x _i

^gleih

y _i

^{sein müssen.}

x _i

^{steht für die} Stützstellen, deren Anzahl

n

^{ist (siehe oben). Damit hat}

r n

Nullstellen. Es gibt nur ein Polynom, dessen Grad kleiner ist als die Anzahl seiner Nullstellen und das ist

die Nullfunktion. Deshalb gilt

r(x) = 0

^{woraus folgt}

p(x) = q(x)

^{und bewiesen}

wäre, dass nur ein Polynom

p(x)

^vomGrad

(n − 1)

^{existiert, welhes durh alle}

Stützstellenverläuft.

DurhInterpolationmitdemPolynomlassensihdieWertezwishendentabella-

risiertenStützstellenapproximieren.Umeinesehr hoheGenauigkeitzu erhalten,

müsste nun eigentlihnur dieAnzahlan Stützstellenentsprehend hohgewählt

werden.Wenn alleStützpunktegleihverteiltwerden ergibtsihjedoheinesehr

groÿeAbweihung an den Ränderndes inder Tabelleangegebenen Intervalls.

DiesesBildverdeutliht dieAbweihung beieinerAnnäherung andieSinuskurve

(gestrihelt) miteinem Polynom bei 20gleihverteiltenStützstellen:

Um diese Abweihung zu beheben ist es sinnvoll die Verteilung der Stützstellen

zu verändern, indem man mehr Stützstellen an den Rändern wählt. Um diese

Verteilung zu optimieren wurde das sogenannte Tshebyshev-Verfahren entwi-

kelt.

BeidiesemVerfahrenwirdnaheinemPolynom

n

^{-tenGradesgesuht,dasindem}

Intervall

[ − 1; 1]

^{denmaximalenBetrag}

1

^{hat sowie alleMaximaund Minimabei}

± 1

^.Dabeisteht

n

^{fürdieAnzahlderspäteren}Stützstellen.DasPolynom,welhes diese Bedingungen erfülltheiÿt Tshebyshev-Polynom

T n

^:

n = 0 : T 0 = 1 n = 1 : T 1 = x n = 2 : T 2 = 2x ² − 1

.

.

.

n = k : T k = cos(k · arccos(x))

Um die Nullstellen zu bestimmen wird

T _n

^{gleih Null gesetzt, da der Cosinus}

seine Nullstellenbei

π

2 + (π)i

^{für alle}

i ∈ N

hat, ergibtsih:

k(arccos(x)) = π

2 + (π)i arccos(x) = π

k ( 1 2 + i) x = cos( π

k ( 1 2 + i))

DerCosinusnimmtaufdemIntervall

[0; π]

^{dieWerteausdemIntervall}

[ − 1; 1]

^an.

DarausergibtsihdieUngleihung

0 ≤ ^π k ( ¹ ₂ +i) ≤ π

^,dasheiÿt

i

^{nimmtalleWerte}

zwishen

i = 0

^und

i = k − 1

^{an. Diesbedeutet auh,dass esingesamt}

k

^Nullstel-

lengibtwelhealsStützstellendienen.DiePolynominterpolationeinerbeliebigen

Funktionandiesen Stützstellenweist eine bedeutendhöhere Genauigkeitauf als

mitgleihverteiltenStützstellen.Eine Vergleihsrehnung zu den gleihverteilten

Stützstellen zeigte nur noh Abweihungen in der Gröÿenordnung der Mashi-

nengenauigkeit und die Graphen der Funktion und der Interpolierenden waren

dekungsgleih.

Wiebereits gesehen, basieren viele Näherungsverfahren auf Polynomauswertun-

gen (Taylorpolynome, Interpolationspolynome). Dabei kommt der Ezienz und

Rehengeshwindigkeit eine immer gröÿere Bedeutung zu. Dies hat zur Folge,

dass versuht wird, dieAlgorithmen zu optimieren.

Das Horner-Shema ist eine spezielle Anordnung der Faktoren eines Polynoms

p(x) = a 0 + xa 1 + x ² a 2 + x ³ a 3 ... + x ⁿ a n

BeiBetrahtungergibtsih,dass hier

n

^{Additionenbenötigtwerden.Hinzukom-}

mennohzurAuswertungder

x ^k

^{-Potenzenjeweils}

(k − 1)

Multiplikationen.Durh densteigendenExponenten addiertmitdenMultiplikationender

x

^-Potenzenund

der dazugehörigen Koezienten

a k

^{ergeben sih insgesamt}

n(n+1)

2

Multiplikatio- nen.

Wirddas

x

ausgeklammert ergibtsih folgedes Polynom:

p(x) = a ₀ + x(a ₁ + xa ₂ + x ² a ₃ ... + x ^{n − 1} a _n )

Nunwird

x

^imentstandenenFaktor wiederausgeklammertsodass sihfolgendes Polynom ergibt:

p(x) = a 0 + x(a 1 + x(a 2 + xa 3 ... + x ^{n − 2} a n ))

DiesesPrinziplässtsih

n

^-malwiederholen,sodasssihfolgendesPolynomergibt:

q(x) = a ₀ + x(a ₁ + x(a ₂ + x(a ₃ ... + xa _n ))...))

DiesesPolynom

q(x)

^{,das durh} äquivalenteUmformungaus dem ursprünglihen Polynom

p(x)

hervorgegangen ist wird jetzt noh auf Ezienz untersuht. Ge- naues Zählen liefert, dass

n

^{Additionen und nur noh}

n

Multiplikationen nötig sind.

6 Fazit

Wozu nun selber rehnen? Mit unserem Projekt sollte verdeutliht werden, wie

ein Computer rehnet und vor allem sollte gezeigt werden, wo die Grenzen un-

serer Rehner liegen. Wenn man diese Grenzen kennt, weiÿ man auh, ab wann

manselbst rehnen muss. Auhwenn dieMashinengenauigkeit beieinem64bit-

Zahlensystem sehr klein ersheint, kann in längeren Rehnungen mitvielen Ite-

rationenso oft aufsummiert werden, dass das Ergebnis entweder verfälsht oder

gänzlihungenauwird.Wenn manjedohweiÿ, dassdasErgebniseinesRehners

einfahwie möglih zu halten, damit dievermutete Ungenauigkeitso gering wie

möglihwird und besser berüksihtigt werden kann.