Re: Stellenauslöschung ! SetPrecision und Algorithmus [Davidson]

["Jens-Peer Kuska" <kuska@XXXXXXX.de>]

Ach so, und wenn man schon Zufallszahlen erzeugt, 
so nützt es wenig,
mit SetPrecision[] diesen Zufallszahlen lauter 
Nullen anzuhängen,
sondern man sollte dann auch die Zufallszahlen mit 
der genwünschten Präzision
erzeugen, dann klappt es auch prima mit 50 
Stellen.

Gruß
 Jens

----- Original Message ----- 
From: "Marc von Bredow" <mvb@XXXXXXX.de>
To: <demug@XXXXXXX.ch>
Sent: Wednesday, May 25, 2005 1:19 PM
Subject: Stellenauslöschung ! SetPrecision und 
Algorithmus [Davidson]


> Hallo,
>
> hier noch einmal die Problematik von N und 
> SetPrecision und die dann
> eintretende Stellenauslöschung beim Rechnen
>
> Das nb beinhalted eine simple Form des 
> Davidson-Verfahren zur Bestimmung des
> betragsmässig größten Eigenwertes einer 100x100 
> Matrix, reell, nicht
> symmetrisch, zufällig.
> Schaut Euch mal die Ergenisse an [mit 
> Mathematica 5.1 gemacht. MIt Version 4
> läuft es so nicht, da ich die Funktion "Norm" 
> und "Eigenvalues" zur
> Bestimmung eines EW benutze. "Eigenvalues[A,1]" 
> liefert den gößten EW. Aber
> das ist im Grunde nebensächlich.
>
> In der ersten Zelle werden die Ausgangsmatrix A 
> und der Startvektor mit N
> berechnet:
>
> n = 100;
> tmp = N[Table[If[i ==  j, Random[Real, {0, 1}], 
> Random[Real, {0, 1}]], {i,
> 1, n}, {j, 1, n}]];
> A = tmp;
> .
> .
> .
> tmp = N[Table[Random[Real, {-10, 10}], {i, 1, 
> n}]];
>
> [Das If ist nebensächlich und dient nur dazu, im 
> Testfall eine
> diagonaldominante Matrix bei Bedarf zu erzeugen, 
> also eine ganz normale
> Matrix mit Reals].
>
> In der zweiten Zelle wird das Gleiche mit 
> SetPrecision mit 50 Stellen
> gemacht
>
> n = 100;
> tmp = SetPrecision[Table[If[i == j, Random[Real, 
> {0, 1}], Random[Real, {0,
> 1}]], {i, 1, n}, {j, 1, n}], 50];
> A = tmp;
> .
> .
> .
> tmp = SetPrecision[Table[Random[Real, {-10, 
> 10}], {i, 1, n}], 50];
>
> Schaut man sich die Ergebnisse an, dann sieht 
> man, daß in der ersten Zelle
> alles gut geht.
> Die MaschinePrecision bleibt während der 
> Rechnung erhalten und die
> Ergebnisse sind in Ordnung.
>
> In der zweiten Zelle wird die Precision 
> aufgefressen, wie man deutlich
> sieht!
> Woran liegt das?
>
> Im Algorithmus wird auch GramSchmidt zur 
> Orthogonalisierung verwendet. Er
> meckert in der zweiten Zelle rum
> \!\(GramSchmidt::"zeromag" \(\(:\)\(\ \)\) "The 
> magnitude of \!\({
>      0``0.51459253908149,
>     0``0.655168185198406, 0``0.6136232110037116, 
> \
> 0.2911151519510895014`0.02911543173540807, 
> 0``0.7102477517539784, \
> 0``0.656413398009062, 0``0.6835108328676073, 
> 0``0.6192498952159269, \
> 0``0.686386657560155, 0``0.7548169826639913, 
> \(\(\[LeftSkeleton] 90 \
> \[RightSkeleton]\)\)}\) was zero under the inner 
> product \!\(Dot\). This \
> suggests that the vectors are linearly 
> dependent."\)
>
> ok.ok, soll uns nicht weiter jucken, das 
> passiert nun mal im Algorithmus.
> Die Vektoren werden tatsächlich nahezu linear 
> abhängig. (Stellenauslöschung)
> Um die Basis neu aufzubauen muß eben genau 
> deswegen der neue Vektor zu den
> bereits bestehenden orthogonalisiert werden. 
> Deswegen ja GramSchmidt, man
> kann auch Householder nehmen, ändert aber nix 
> wesentliches.
> In der ersten Zelle klappt's, in der zweiten 
> nicht.
> Warum?
>
> Danke und beste Grüße,
> Marc von Bredow