| Autor |
 cos(x) mit über 10 Mio. Stellen |
|
hyperG
Senior 
Dabei seit: 03.02.2017
Mitteilungen: 2159
 |  Themenstart: 2017-11-22
|
Natürlich könnte man teure Software wie Mathematica verwenden,
ABER ich möchte das z.B. mit Pari/GP.
gegeben:
x=0.9528... (über 10 Mio. weitere Stellen)
gesucht: cos(x)
Problem: aktuelle Rechenzeit über 16 Stunden & nachts muss der PC aus sein.
Pari/GP kann parallel nur mit Vektoren (Arrays), aber hier liegt als Quelle 1 reelle Zahl in einer Datei.
Die Gesamtzeit muss nicht sehr viel schneller werden, sondern es reicht, wenn die größte Teilzeit unter 10 h kommt.
Grundrechenarten sind schnell genug.
Idee Nr 1:
Wegen der Reihenentwicklung und Unabhängigkeit von Pi sollte es mit Winkelhalbierung schneller gehen:
cos(x)=1-2*sin(x/2)²
dann:
sin(5x)=sin(x)*(5-20*sin(x)²+16*sin(x)^4)
cos(x)=1-2[sin(x/10)*(5-20*sin(x/10)²+16*sin(x/10)^4)]²
Ich muss also nur 1 mal sin(x/10) berechnen (1 Kommastelle verschieben), was schon mal schneller werden sollte.
Andere Ideen?
Idee Nr. 2:
hypergeometrische Funktionen
cos(x)=hyg0F1(1/2, -x²/4)
Weitere Ideen??
|
 Profil
|
Primentus
Senior 
Dabei seit: 18.02.2016
Mitteilungen: 2180
Wohnort: Deutschland
 | Beitrag No.1, eingetragen 2017-11-22
|
Hallo hyperG,
einen Algorithmus, wie man möglichst schnell einen 10millionenstelligen Cosinuswert aus einem 10millionenstelligen Argument berechnet, hab ich zwar nicht parat, möchte aber trotzdem mal einwerfen, dass Mathematica für diese Berechnung lediglich 2 Minuten benötigt.
Wie wird denn Dein Argument 0.9528... genau erzeugt? Liegt das schon vor oder ist das auch Teil des Berechnungsprozesses? Ich vermute mal, dass bei Deiner bisherigen Berechnung die Erzeugung dieses Argumentes am meisten Zeit beansprucht.
LG Primentus
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.2, vom Themenstarter, eingetragen 2017-11-23
|
Hallo Primentus,
ich hätte statt "über" besser "weit über" schreiben sollen.
Die fast 17 h bezogen sich auf 100 Mio. Stellen.
Bei 2 Mio. Stellen schafft GP/Pari mit 1 Kern 97 s ... 102 s (i9...i7).
Mathematica benutzt ja schon in der 350 € Version 4 Kerne.
Der Anstieg der Rechenzeit mit der Stellenanzahl ist stark nichtlinear.
Ein Versuch mit Verkleinerung des Argumentes (Idee 1 ;selbst mit Faktor oberhalb 30000) brachte bei 100 s nur 0,5 s für die 1. Teilrechnung. Die Gesamtrechenzeit änderte sich kaum.
Vermutlich hat Pari bereits viel optimiert...
Ich stehe nun vor der Entscheidung, mir entweder selbst was mit "parallel" zu schreiben, oder was teures zu kaufen...
Oder jemand hat noch andere Ideen...
Oder ich zerlege selbst die sin-Funktion in eine Reihe,
und berechne jeden Tag ein Teilbereich der Reihe.
Am letzten Tag werden alle Teilsummen (die ja in Dateien vorliegen), zusammenaddiert.
Das müsste sogar parallel funktionieren, da man GP/Pari parallel starten kann...
|
Profil
|
mrdjv2
Aktiv 
Dabei seit: 05.07.2003
Mitteilungen: 1020
Wohnort: Aachen
 | Beitrag No.3, eingetragen 2017-11-23
|
Hallo,
ich hab zwar von Numerik und so nur relativ wenig Ahnung, aber ich wüsste schon gern, was du genau möchtest.
Hast du eine bestimmte Anzahl von Nachkommastellen, die du berechnet haben möchtest? Wenn ja, wie viele?
"Weit über 10 Mio" ist immer so eine Sache...wenn du schon die 100 Mio "geknackt" hast, geht es dir dann darum, die 100 Mio schneller zu bekommen oder möchtest du noch mehr Nachkommastellen haben?
Viele Grüße
mrdjv2
|
Profil
|
MontyPythagoras
Senior 
Dabei seit: 13.05.2014
Mitteilungen: 3417
Wohnort: Werne
 | Beitrag No.4, eingetragen 2017-11-23
|
Hallo hyperG,
bringt vielleicht das unendliche Produkt etwas?
\[\cos x=\prod_{n=1}^{\infty}\left(1-\frac{4x^2}{\pi^2\left(2n-1\right)^2}\right)\]
Ciao,
Thomas
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.5, vom Themenstarter, eingetragen 2017-11-23
|
Hallo mrdjv2,
Problem ist die Zeit, wie lange mein PC eingeschaltet werden darf:
13 Stunden
Rechenleistung habe ich genug: 20 + 8 CPU-Kerne
Da 99% der Software immer nur 1 Kern ansprechen, werde ich innerhalb der 13 h nicht fertig.
Lösung 1:
Summe, denn die kann man gut parallelisieren und in beliebig kleine Teilsummen aufteilen.
Habe bereits Argument bis 0.0000165767... (und 100 Mio. weitere Stellen)
verkleinert.
Bei der gewünschten Genauigkeit macht das 7.72 Mio. Summanden
könnte man auf 20 EXE mit je 386000 Summanden aufsplitten
Hallo Thomas,
Das Partialprodukt Deines Vorschlages lautet
Pi*Pochhammer(1/2 - x/Pi,n)Pochhammer(x/Pi + 1/2,n)/Gamma(n + 1/2)^2
Das x² sieht zunächst gut aus, aber 2 große Nachteile:
\sourceon Tabelle der Genauigkeit:
Konvergenzgeschwindigkeit (Größe von n):
n Genauigkeit bei x=1/100 x=1/1000 0.000006098227
17 -5.9e-7
1000 -1.01e-8 -1e-10
1e6 -5e-11 -1.5e-13 -1.9e14
\sourceoff
also bei 1 Mio. Faktoren nur 14 Stellen Genau :-(
Und Fehlerfortpfanzung bei Mul ist langsamer & größer als Add...
|
Profil
|
mrdjv2
Aktiv
Dabei seit: 05.07.2003
Mitteilungen: 1020
Wohnort: Aachen
| Beitrag No.6, eingetragen 2017-11-23
|
Ich weiß immer noch nicht, was du willst.
Möglichst viele Nachkommastellen in einer bestimmten Zeit oder eine feste Anzahl an Nachkommastellen in möglichst kurzer Zeit??
|
Profil
|
Primentus
Senior
Dabei seit: 18.02.2016
Mitteilungen: 2180
Wohnort: Deutschland
| Beitrag No.7, eingetragen 2017-11-23
|
@hyperG:
Lange Zeit habe ich ja einen inzwischen schon älteren Rechner mit 32 Bit-Betriebssystem und einer ziemlich alten Mathematica-Version genutzt, die für heutige Verhältnisse ebenfalls recht langsam war. In letzter Zeit habe ich aber nochmal einiges investiert und mir einen schnelleren Rechner mit immerhin 4 Kernen sowie 64 Bit-Betriebssystem und einer aktuellen Mathematica-Version zugelegt. Sicherlich überlegt man es sich gut, ob man diese knapp 370 Euro allein für Mathematica ausgeben will, aber ich kann nur sagen, dass dieses sehr leistungsfähige Programm absolut jeden Euro Wert ist. Mit der neuen Version sind die Berechnungszeiten bei einer Vielzahl von Berechnungen um ein Vielfaches besser als bei der alten Version. Die Tatsache, dass dort 4 Kerne arbeiten, macht sich auf jeden Fall bemerkbar. Und man kann hier wirklich mit recht wenigen Codezeilen zu umfangreichen Ergebnissen kommen. Sicherlich kann man aber auch ein anderes CAS verwenden. Ich persönlich möchte nur sagen, dass ich mit Mathematica durchwegs sehr gute Erfahrungen gemacht habe.
Die von mir angegebenen 2min Rechenzeit bezogen auf 10 Millionen Stellen beziehen sich bereits auf mein neues System. Habe jetzt, da Du von nicht-linearem Anstieg gesprochen hast, mal für ein 100millionenstelliges Argument den 100millionenstelligen Cosinuswert ermittelt. Dabei ergab sich eine Rechenzeit von ca. 41,5 Minuten (vorausgesetzt, dass das zu übergebende Argument bereits fertig vorliegt). Nicht-linearer Anstieg stimmt schon mal, da es bei linearem Anstieg um die 20 Minuten hätten sein müssen. Aber ich schätze mal, dass man innerhalb Deiner 13 Stunden-Vorgabe mit Mathematica schon mehrere hundert Millionen Stellen berechnen könnte. Vielleicht ist es also doch eine Überlegung wert, dass Du Dir ein solches kommerzielles System zulegst.
Trotzdem bin ich aber auch noch gespannt, ob jemand noch eine Alternativberechnung in petto hat, mit der es noch schneller geht.
LG Primentus
|
Profil
|
Delastelle
Senior
Dabei seit: 17.11.2006
Mitteilungen: 2515
 | Beitrag No.8, eingetragen 2017-11-23
|
Hallo hyperG!
Die 16 Stunden vom Themenstart können nicht in 2x 8 Stunden mit Zwischenspeicherung und wieder aufsetzen auf dem Zwischenergebnis umgewandelt werden?
Viele Grüße
Ronald
|
Profil
|
rlk
Senior
Dabei seit: 16.03.2007
Mitteilungen: 11660
Wohnort: Wien
 | Beitrag No.9, eingetragen 2017-11-24
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.10, vom Themenstarter, eingetragen 2017-11-24
|
Hallo Ronald,
ja, solche Algorithmen mit quadratischer Konvergenzgeschwindigkeit sind für große Zahlen geeignet.
Allerdings habe ich noch keine konkrete Iteration (nur Näherungen für unter 100 Stellen) gefunden, die ich hier anwenden könnte.
@Delastelle:
Das ist es ja gerade: wenn man Fremdsoftware einsetzt, muss man entweder abbrechen (alles umsonst), oder bis zum bitteren Ende warten.
Deshalb hatte ich ja die Idee mit der Reihe, die man beliebig abbrechen kann. Bedeutet aber über 7 Mio. mal je 4 Grundrechenarten mit je über 100 Mio. Stellen selbst verwalten...
Ich habe selbst mit Win10 den "Ruhestand" probiert. Es sah am nächsten Tag auch gut aus, Aber statt 100 Mio. wurde nur 1 Mio. abgelegt???
(Ruhezustand ist was für Editoren, aber nicht für Programme, die viele Init-Parameter haben und 100% pro Kern ziehen.)
@mrdjv2:
ich will ein Ergebnis mit den hier schon 3 mal genannten Randbedingungen, ohne teure Software zu kaufen.
Da Du immer wieder nachfragst eine Analogie:
Du hast eine Supersportwagen (PC), der 350 km/h kann (genug Kerne und 512 Bit Befehle hat)...
darfst aber nur im 1. Gang (1 Kern und vorgegebene Software) fahren (also nur 80 km/h)...
der Weg ist Berg & Tal, aber Du hast keine Bremse und darfst nur tags fahren. Da Du mit 80 km/h nicht den Berg bis zur Nacht schaffst (22:30 Uhr wird PC ausgeschaltet), rollst Du immer wieder zurück. (kein Ergebnis)
Natürlich könntest Du ein teures Elektroauto (Mathematica) kaufen...
Ich versuche nun, höhere Gänge (andere Algorithmen) zu finden. Natürlich ist schneller = besser, ABER auch ein sicherer HALT (Bremse oder Zwischensicherung) führt am nächsten Tag auch zum Ziel.
Erste Erfolge bis jetzt:
Habe im Internet eine 64Bit DLL gefunden, die den Befehl
mpfr_sin() enthält und bekam das STAUNEN:
Bei 670000 Stellen (!), wo GP/Pari 16,58 s braucht,
stand das Ergebnis schon nach 0,843 s fest!!!
Fast 20 mal schneller -> das musste überprüft werden:
- Übereinstimmung nur bis Stelle 201702
- noch größere Werte führten leider zum Absturz??? (RAM habe ich genug)
Schade :-(
Aber so leicht gebe ich nicht auf. Bei vielen meiner Projekte war ich am Ende immer schneller als Mathematica...
|
Profil
|
AlphaSigma
Senior
Dabei seit: 23.11.2012
Mitteilungen: 462
| Beitrag No.11, eingetragen 2017-11-24
|
Hallo hyperG,
nur aus Interesse: Was ist der Grund für die Randbedingung
"nachts muss der PC aus sein" und "Problem ist die Zeit,
wie lange mein PC eingeschaltet werden darf: 13 Stunden"?
Ist das eine künstliche Randbedingung, um eine Herausforderung zu haben
oder steht der Rechner im Schlafzimmer?
|
Profil
|
rlk
Senior
Dabei seit: 16.03.2007
Mitteilungen: 11660
Wohnort: Wien
| Beitrag No.12, eingetragen 2017-11-25
|
Hallo hyperG,
Algorithmen zur Berechnung von Winkelfunktionen werden in dem Artikel
Fast Multiple-Precision Evaluation of Elementary Functions
von Richard P. Brent angedeutet. Du findest ihn zum Beispiel hier:
http://maths-people.anu.edu.au/~brent/pub/pub034.html
http://cr.yp.to/bib/1976/brent-elementary.pdf
Kannst Du etwas mehr zum Hintegrund der Berechnung schreiben?
Servus,
Roland
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.13, vom Themenstarter, eingetragen 2017-11-26
|
Hallo Roland,
die atan-Iteration habe ich hinbekommen und konnte sie per Iterationsrechner ausprobieren (Beispiel 139):
Iterationsrechner
\sourceon Ergebnistabelle des Iterationsrechners
i aB aC aD
0 0.4001419795766517 1.1920928955078122e-7 -4.440892098500626e-16
1 0.19264589033659404 0.0006905338836844341
2 0.09551017770267872 0.05251975614959004
3 0.05013807178516839 0.4354728560800887
4 0.03595650152019444 0.9194234067524564
5 0.03450449451677201 0.9991184687781045
6 0.03448924990605769 0.9999999027771472
7 0.03448924822549272 0.9999999999999989
\sourceoff
in aD[0] steckt die Differenz zu atan(0.952847864654941945)
ABER: in der Abschlussberechnung steckt log(x) mit drin,
die leider - wie alle trigonometrischen Funktionen -
nicht über die 200000 Stellen kommt.
Nur Grundrechenarten und Wurzel sind schnell & kommen über 100Mio.
In der PDF steht auch was von
T(m)=tan(...???
aber was ich auch probierte, es kam nichts bekanntes heraus.
Nur bei m=0.5 {also im Iterationsrechner a=0.5}
T(m) Iteration
Kam e^(Pi/2) raus {und auch ohne log-Funktion}
Die Wandlung sin(x)=2*tan(x/2)/(1+tan(x/2))
ist ja kein Problem. Aber wie bekomme ich die Iteration hin.
Auch AGM(x,y) {arith.-geom- Mittel ist kein Problem für über 100 Mio. Stellen ist kein Problem}
aber es gibt in der PDF keine Parameter dafür??
Habe auch schon an Iteration gedacht:
x:= x-(asin(x)-a)*sqrt(1-x²)
Aber das Problem ist asin...
Nun könnte man asin per Näherung oder Iteration...
aber man kommt von 100ten zum 1000ten...
Hintergrund sind viele Iterationsberechnungen, die ich oft brauche (Ziegenfaktor usw.)
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.14, vom Themenstarter, eingetragen 2017-11-30
|
Nach so viel Arbeit nur Rückschläge:
- hypergeometrische Funktionen: zu langsame Konvergenz :-(
- Newton-Iteration per asin -> asin per atan Iteration, ABER am Ende wieder log(x) mit der 200000 Begrenzung :-(
- mit den superschnellen Funktionen mpfr_:
als AGM(x,y) {100Mio in 335 s!!} und sqrt() bis 100Mio {Argumente um 1 herum} funktionierte, war ich noch optimistisch;
beim Einlesen von 300000stelligen Dezimalzahlen (Datei-> String; egal wo Komma) war dann Schluss -> selbst mit Quellcode-Umschreibung stieß ich immer wieder auf die Begrenzung von 200000 Stellen! Als selbst die Grundrechenart Mul(x,y) nur mit y vom Typ uint64 funktionierte war's dann endgültig aus. Sobald y vom Typ reelle Zahl verwendet wird, war da wieder die 200000 Stellen Begrenzung. Selbst einfache Rundung bezog sich immer auf Bit-Stelle und nicht auf Dezimalstelle...
So viele Funktionen umgeschrieben und immer wieder neue Begrenzungen
:-(
|
Profil
|
digerdiga
Ehemals Aktiv 
Dabei seit: 15.11.2006
Mitteilungen: 1389
| Beitrag No.15, eingetragen 2017-11-30
|
Vielleicht bin ich gerade zu naiv und versteh das Problem nicht, aber wenn ich in Julia cos(0.95......) berechne braucht der 41 Sekunden für 10^7 Stellen. 100Millionen hab ich jetzt nicht probiert, aber ich denke es geht auch in die Richtung 20-30 Minuten.
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.16, vom Themenstarter, eingetragen 2017-12-02
|
Hallo digerdiga,
danke für den Tipp!
2 sehr interessante Fakten:
1. neue Software analog Pari/GP die auch über 10 Mio. Stellen importieren kann
2. die darin enthaltene DLL ist genau der reparierte Nachfolger meiner fehlerhaften DLL (mpfr_ Funktionen mit Abstürze oberhalb 200000 Dezimalstellen)
Ich kann durch einfaches Austauschen der DLL alle meine selbst optimierten Funktionen auch oberhalb 10 Mio. Stellen nutzen! Arbeit war nicht umsonst.
Hinweis: precision bezieht sich bei mpfr_ Funktionen auf Bit-Genauigkeit. Um "echte Dezimalstellen" zu bekommen, muss Faktor log(10)/log(2) hinzugefügt werden, was dann statt der 41 s (10^7 Stellen) bei mir etwa 97 s ausmacht. Außerdem ist Argument nicht 0.95, sondern größer und kommt aus einer Datei mit 100 Mio Stellen.
ABER nun zur Überprüfung der Ergebnisse:
\sourceon cos(x/2)*2 10 Mio. Stellen
Software Rechenzeit richtige dez. Stellen
Pari/GP 1433 s 10 Mio.
Mathematica 120 s 10 Mio.
Julia 97 s 306916
cpp 73 s 315096
\sourceoff
Bei Pari ist es egal, ob cos(x/2... oder sin(x/312500 gerechnet wird, da die 1 s Zeiteinsparung bei sin bei der nachträglichen Rückrechnung wieder aufgebraucht wird.
Da ich den Code der mpfr_ Funktionen nun halbwegs kenne, könnte es 4 mögliche Ursachen für die "nur 300000" richtige Stellen {trotz setprecision(33219272)} geben:
- strInput kann Genauigkeit nicht auf die Variable übertragen
- intern gibt es "GetPi()" (und log2), was begrenzt sein kann
- der Algorithmus bricht zu früh ab
- Fehlerfortpflanzung
Bevor die Ergebnisse bei 10 Mio. nicht stimmen, brauche ich nicht die 100 Mio. starten...
Aber ich kann ja nun mit meiner Newton-Iteration {beinhaltet asin}
noch einiges optimieren...
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.17, vom Themenstarter, eingetragen 2017-12-02
|
Ich hab's !!!!!
statt
mpfr_sin oder _cos oder
mpfr_sincos_fast
bin ich noch eine Quellcode-Ebene herunter und bei Anpassung von
sincos_aux {ohne Hilfskonstanten + Array-Vergrößerung}
stimmten nun ALLE Dezimalstellen!!!
Nun also 10 mal mehr -> Start mit 100 Mio:
[precision is 332193078 bits]
time = 1228.595 s also unter 21 min !!!
48 mal schneller als Pari/GP !!!
2 mal schneller als 4-Kern-Mathematica!
Bis 50 Mio stimmen alle Dezimalstellen schon mal.
Nun warte ich auf Primentus, was Mathematica sagt...
und das mit i7 & 1 Kern -> mit i9 & 20 Kernen sind da noch gewaltig Optimierungsmöglichkeiten...
Zwischenergebnisse kann ich jederzeit in Dateien speichern!
Und wenn man andere Algorithmen ansetzt (dabei nur die funktionierenden mpfr_ Funktionen verwendet), kann man auch selbst überprüfen, bis zu welcher Dezimalstelle alles stimmt.
Bin gespannt, wann nächste Obergrenzen (Überläufe, Abstürze) kommen.
Danke an alle, die geholfen haben!
Grüße Gerd
|
Profil
|
digerdiga
Ehemals Aktiv 
Dabei seit: 15.11.2006
Mitteilungen: 1389
| Beitrag No.18, eingetragen 2017-12-02
|
\quoteon(2017-12-02 12:57 - hyperG in Beitrag No. 17)
Ich hab's !!!!!
statt
mpfr_sin oder _cos oder
mpfr_sincos_fast
bin ich noch eine Quellcode-Ebene herunter und bei Anpassung von
sincos_aux {ohne Hilfskonstanten + Array-Vergrößerung}
stimmten nun ALLE Dezimalstellen!!!
Nun also 10 mal mehr -> Start mit 100 Mio:
[precision is 332193078 bits]
time = 1228.595 s also unter 21 min !!!
48 mal schneller als Pari/GP !!!
2 mal schneller als 4-Kern-Mathematica!
Bis 50 Mio stimmen alle Dezimalstellen schon mal.
Nun warte ich auf Primentus, was Mathematica sagt...
und das mit i7 & 1 Kern -> mit i9 & 20 Kernen sind da noch gewaltig Optimierungsmöglichkeiten...
Zwischenergebnisse kann ich jederzeit in Dateien speichern!
Und wenn man andere Algorithmen ansetzt (dabei nur die funktionierenden mpfr_ Funktionen verwendet), kann man auch selbst überprüfen, bis zu welcher Dezimalstelle alles stimmt.
Bin gespannt, wann nächste Obergrenzen (Überläufe, Abstürze) kommen.
Danke an alle, die geholfen haben!
Grüße Gerd
\quoteoff
Womit hast du das jetzt gelöst? c++?
Doofe Frage aber wie hast du genau überprüft, dass z.B. bei Julia die Stellen nach 300,000 falsch sind?
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.19, vom Themenstarter, eingetragen 2017-12-02
|
\quoteon(2017-12-02 13:43 - digerdiga in Beitrag No. 18)
...
Womit hast du das jetzt gelöst? c++?
\quoteoff
Ja, c++ Visual Studio Express, weil der Quellcode aller mpfg_ Funktionen auch in c geschrieben wurde. Anstatt jedoch Mio. Zeilen von Code (der für LINUX geschrieben wurde) neu für Win lauffähig zu machen (und hunderte Compilerschalter anzupassen), klinke ich mich in die fertige DLL ein -> genau so wie Julia das auch macht.
\quoteon(2017-12-02 13:43 - digerdiga in Beitrag No. 18)
... wie hast du genau überprüft, dass z.B. bei Julia die Stellen nach 300,000 falsch sind?
\quoteoff
Wie schon angedeutet gibt es mehrere Wege für Validierungsberechnungen:
- Vergleich mit anderer Software (Pari/GP, bc , Mathamatica,...)
- Vergleich mit anderen Funktionen:
cos(x)=1-2*sin(x/2)²
cos(x)=1-2[sin(x/10)*(5-20*sin(x/10)²+16*sin(x/10)^4)]²
- Vergleich mit anderen Algorithmen:
cos(x)=hyg0F1(1/2, -x²/4)
Selbstkonvergenz:
x:= x-(asin(x)-a)*sqrt(1-x²)
Start mit x[0]=sin(..) aber nur 200000 Stellen (bis dahin noch OK)
danach kann man asin mit der atan-Funktion berechnen (für die es auch Iterationen gibt).
Die Newton-Iterationen konvergieren schnell: Verdopplung der richtigen Stellenzahl pro Schritt -> 400000 Stellen usw.
Das Komplizierte bei den mpfr_ Funktionen ist nun, dass alles auf Binär optimiert ist. Man muss Ergebnisse (oder Zwischenergebnisse) in Strings (txt Datei) ablegen.
Die neue DLL stürzt zwar nicht mehr ab, aber man muss jede Funktion überprüfen: wenn dann der Byteweise Vergleich Unterschiede 2er Dateien
(Referenz muss nicht mal zu 100% stimmen)
zu Tage bringt, kann was nicht stimmen!
Bei den Iterationen muss man unterscheiden zwischen:
- festen exakten Startwerten z.B. agm(1,x)
{hier bewirkt eine kleine Störung sofort falsches Ergebnis)
- selbstkonvergierende Iteration
{hier kann man erste Schritte mit weniger Stellen rechnen und erst in letzter Iteration voll hochschalten -> wie ein Regler pendelt sich die Iteration auf höchste Genauigkeit ein}
Beispiel: als ich Julia & meine cpp Optimierung verglichen habe, gab es schon ab Position 306916 Unterschiede. Ein weiterer Vergleich mit Pari/GP zeigte bei dem cpp Ergebnis ab Pos 315096 Unterschiede.
Die Pari-Ergebnisse hatte ich per Selbstkonvergenz bereits bis 50 Mio validiert.
Dann gibt es noch BBP-Formeln, wo man einzelne Hex-Bytes vergleicht. Der Weltmeister der Pi-Software überprüft so etwa alle 1 Mrd. Stellen, ob noch alles stimmt. Leider sind noch keine BBP-Formeln für Dezimale Stellen bekannt. Wenn also bei der Wandlung von Hex nach Dez Fehler passieren, ist Ergebnis Falsch und wurde nicht für Dezimalstellen validiert.
Da die letzte Berechnung bis 50 Mio. validiert ist, warte ich nun auf Primentus, dem ich meine letzten 300 Stellen schickte.
Wenn ich mehr Zeit hätte, könnte ich auch selbst per oben beschriebenen Wege weitere Berechnungen starten (die das selbe Ergebnis bringen sollten) und neue Byte-Vergleiche starten.
Grundrechenarten, Wurzel & AGM scheinen die einzigen mpfr_ Funktionen zu sein, die um 1 herum keine Begrenzung im Bereich 150000...320000
haben. Weicht man stark vom Wert 1 ab (was bei Zwischenergebnissen schnell passieren kann) verschiebt sich nochmal alles. (so wie die Wandlung in Strings, die ich Blockweise mit 200000 Stellen realisieren wollte) Da das Runden auch auf Binärzahlen bezogen war und ich den Umweg über die Multiplikation mit sehr großen Zahlen gehen wollte, gab es auch Abstürze.
Bei den Mio. Zeilen hardwarenaher c Sprache ist Fehlersuche sehr zeitaufwendig!
Bei einer FFT-Multiplikation hatte ich mal 3 Wege für die sin- Zwischenergebnisse: FPU, CPU, AVX-Befehle -> heraus kamen 3 unterschiedliche Maximalgenauigkeiten & Geschwindigkeiten (winzige Unterschiede in der versteckten 16. Dezimalstelle, da binär).
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.20, vom Themenstarter, eingetragen 2017-12-03
|
Einen 2. Kern konnte ich an einigen Stellen noch aktivieren.
So wurde aus
\quoteon(2017-12-02 00:27 - hyperG in Beitrag No. 16)
\sourceon cos(x/2)*2 10 Mio. Stellen
Software Rechenzeit richtige dez. Stellen
Pari/GP 1433 s 10 Mio.
Mathematica 120 s 10 Mio.
Julia 97 s 306916
cpp 73 s 315096
\sourceoff
\quoteoff
Neu:
\sourceon cos(x/2)*2 10 Mio. Stellen
Software Rechenzeit richtige dez. Stellen
Pari/GP 1433 s 10 Mio.
Mathematica 120 s 10 Mio.
cpp 2 threads 62 s 10 Mio.
\sourceoff
(beim i9 nur 60 s)
Und meine Ausgangsrechnung mit 100 Mio dezimalen Nachkommastellen:
\sourceon cos(x/2)*2 100 Mio. Stellen
Software Rechenzeit richtige dez. Stellen
Pari/GP 16,8 h 100 Mio.
Mathematica 41,5 min 100 Mio.
cpp 2 threads 16,8 min 100 Mio.
\sourceoff
Also Faktor 60 -> na wenn das kein Erfolg ist :-)
|
Profil
|
digerdiga
Ehemals Aktiv
Dabei seit: 15.11.2006
Mitteilungen: 1389
| Beitrag No.21, eingetragen 2017-12-03
|
Ich möchte das mal selber ausprobieren:
Zu Testzwecken reicht es ja wenn man cos(0.9528) bis auf 10^7 Stellen berechnet denke ich. Was sind deine letzten 100 Stellen in diesem Fall?
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.22, vom Themenstarter, eingetragen 2017-12-03
|
Ergebnisse:
\sourceon cos(0.9528)
# Soll ab Pos 306890 459423813195289597370437889596885917728725...
# Julia Pos 306890 459423813195289597370437889596885962619893...
# Soll Pos 9999980 875818151444493976476...
# Julia Pos 9999980 20742729868253737555...
\sourceoff
ab Position 306924 alles falsch! (nicht mal 5% richtig)
Auch
asin(x) sah bis 10 Mio Stellen noch gut aus, bei 100 Mio Absturz...
Validierung ist oberhalb 200000 Stellen absolut wichtig!
|
Profil
|
digerdiga
Ehemals Aktiv
Dabei seit: 15.11.2006
Mitteilungen: 1389
| Beitrag No.23, eingetragen 2017-12-04
|
Das bekomm ich auch (bei den 10^7 allerdings ganz zum Schluss: ...7375566) . Was ich dabei nicht ganz verstehe ist was der da macht!?! Man würde doch erwarten, dass die das richtig implementiert haben...
|
Profil
|
Primentus
Senior
Dabei seit: 18.02.2016
Mitteilungen: 2180
Wohnort: Deutschland
| Beitrag No.24, eingetragen 2017-12-04
|
\quoteon(2017-12-03 20:50 - hyperG in Beitrag No. 22)
Ergebnisse:
\sourceon cos(0.9528)
# Soll ab Pos 306890 459423813195289597370437889596885917728725...
# Julia Pos 306890 459423813195289597370437889596885962619893...
# Soll Pos 9999980 875818151444493976476...
# Julia Pos 9999980 20742729868253737555...
\sourceoff
\quoteoff
Nur zur Sicherheit als Info:
Die beiden Soll-Werte stimmen exakt mit den Ziffernfolgen überein, die auch Mathematica liefert.
LG Primentus
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.25, vom Themenstarter, eingetragen 2017-12-04
|
Validierungsberechnungen abgeschlossen:
sin_cos_aux
stimmt nun mit der selbstkonvergierenden Iteration
(siehe Beitrag 13 )
also
(x-(acos(x)-1.9056957293.../2)*sqrt(1-x²))*2,x=1.15872847301812.../2
auf 100 Mio. + 80 Stellen überein!
Allerdings musste ich acos(x)=atan(sqrt(1-x²)/x)
ersetzen, da mpfr_acos zw. 10 ... 100 Mio abstürzte. (vermutlich wieder die intern gespeicherte Konstante Pi oder Array-Grenzen)
Trotz des Umweges lag das Ergebnis in weniger als 31 min vor, was immer noch schneller als Mathematica ist.
Danke Primentus für Deine Bestätigung:
a) Deine Ergebnisse stimmen mit meinen auf 100000080 Stellen!
b) danke auch für die Bestätigung der Soll-Werte
Also eine Doppelt-starke Validierung!
Interessant ist die log(x) Funktion:
\sourceon log(1.158) auf 10 Mio. Stellen
mpfr_log 95.290 s warum so langsam?
Pari/GP 45.599 s
\sourceoff
mpfr_log rechnet mit der Pi/(4AGM(...) Iteration, die vermutlich durch das Pi (was mit berechnet werden muss) so langsam wird...
Bei Pari/GP habe ich den Algorithmus noch nicht gefunden.
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.26, vom Themenstarter, eingetragen 2017-12-04
|
\quoteon(2017-12-04 00:18 - digerdiga in Beitrag No. 23)
... Man würde doch erwarten, dass die das richtig implementiert haben...
\quoteoff
Das war bei der alten DLL ja noch schlimmer: Abstürze
Das bedeutet, dass niemand je über 300000 Stellen damit berechnete!
Absturz-Fehler sind zwar nun behoben, aber wie ich im Beitrag 16 zeigte:
- strInput kann Genauigkeit nicht auf die Variable übertragen
- intern gibt es "GetPi()" (und log2), was begrenzt sein kann
- der Algorithmus bricht zu früh ab
- Fehlerfortpflanzung
strInput & Fehlerfortpflanzung habe ich nun ausgeschlossen.
GetPi(... konnte ich durch Einschränkung der Argumente und Nachprogrammierung der sincos_aux-Funktion umgehen.
Das nächste sind Arrays: intern wird oft mit mpz_ Funktionen gearbeitet -> also Ganzzahlen! Bis zu einer oberen Grenze geht das auch schneller -> aber keiner hat die obere Grenze ausgetestet!
Arraygrenzen sind oft an die 64 Bit gebunden -> die habe ich bei meinen Funktionen vergrößert...
Das ist es ja gerade: nichtlineare Funktionen oberhalb 300000 Stellen haben es in sich!
Selbst meine Funktionen rechnen nur in einem engen Bereich um 1 herum richtig.
Iterationsalgorithmen konvergieren nicht bei allen Argumenten...
Sonst könnte ja jeder täglich neue Rekorde brechen...
hier
|
Profil
|
digerdiga
Ehemals Aktiv
Dabei seit: 15.11.2006
Mitteilungen: 1389
| Beitrag No.27, eingetragen 2017-12-04
|
Welche Version von Julia hast du benutzt?
Scheint wohl an MPFR zu liegen was erst mit 3.1.5 in v0.6.1 implementiert wurde...
https://github.com/JuliaLang/julia/issues/22758
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.28, vom Themenstarter, eingetragen 2017-12-04
|
Die neuste Julia-Version, die man auf der Internetseite bekommt.
Die darin enthaltene DLL ist auch neu: 25.10.2017
Ich glaube nicht, dass in naher Zukunft sich was bei Nachkommastellen oberhalb 300000 ändern wird -> ist bisher (seit Jahr 2000) auch nicht aufgefallen...
Was sich immer mal ändert, sind Geschwindigkeitsoptimierungen oder grobe Fehler.
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.29, vom Themenstarter, eingetragen 2017-12-07
|
Und da sind sie schon die neuen Obergrenzen:
a) bis 200 Mio. noch alles OK: sincos (beide zusammen) nicht mal 38 min
b) ab 250 Mio. schlagen die mpz_ Funktionen zu:
(UNLIKELY (new_alloc > INT_MAX)) ...
Selbst wenn ich das ausbügeln würde (was sehr viel Zeit verschlingt)
kommt schon die nächste Grenze bei einfacher Multiplikation:
c)
mpfr_mul(): macht nicht mehr als
2^(32-2)*log(2)/log(10) = 323228497 (etwas über 300 Mio. dez. Stellen)
Vom Weltmeister könnte ich noch die FFT-Mul. einbauen,
aber die geht auch nur bis etwa 800 Mio. (die sin-Werte in double Genauigkeit geben nicht mehr her)
Natürlich könnte man in Verdopplungs-Schritten die
Karazuba multiplikation verwenden... -> alles viel Arbeit...
|
Profil
|
Carmageddon
Senior
Dabei seit: 22.12.2009
Mitteilungen: 664
| Beitrag No.30, eingetragen 2017-12-07
|
Ich muss einfach mal fragen - Auch im Bezug zu deinen anderen Themen: Warum und für was benötigt man den bitte 200 Millionen Stellen Genauigkeit? Außer zu einem rein akademischen Zweck (reiner Selbstzweck) fällt mir da leider absolut nichts ein? :-?
Sry für Offtopic...
|
Profil
|
digerdiga
Ehemals Aktiv
Dabei seit: 15.11.2006
Mitteilungen: 1389
| Beitrag No.31, eingetragen 2017-12-07
|
Im Rahmen politischer Korrektheit kann man gar nicht exakt genug sein...
Das gilt besonders für die Herren Volksvertreter!
Da ist mir lieber jemand der auf diese Art "exakt" ist!
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.32, vom Themenstarter, eingetragen 2017-12-07
|
\quoteon(2017-12-07 13:34 - Carmageddon in Beitrag No. 30)
... Warum und für was benötigt man den bitte 200 Millionen Stellen Genauigkeit? ... (reiner Selbstzweck) fällt mir da leider absolut nichts ein? :-?
\quoteoff
Normalerweise rechtfertige ich mich nicht für jeden meiner Rekorde:
133731Ziegenfaktor
Ich muss da mal eine Seite aller Begründungen anlegen, damit ich nicht in jedem Forum oder anderen Treffen immer wieder Fakten neu sammeln muss...
Natürlich provozieren Sätze wie: "den kann man nur näherungsweise berechnen..."
Auch habe ich bis heute mein 8stelliges Geburtsdatum darin noch nicht gefunden...
Bei Pi waren es nicht mal 70 Mio. siehe
Pi Datenbank
Interessant sind Vergleiche mit anderen Ziffernhaeufigkeiten wie Vergleich Ziffernhaeufigkeit mit Pi
... und natürlich die unterschiedlichen Algorithmen, die alle zum selben Ergebnis kommen. Pi ist "ausgelutscht", da kann man als Privatperson nichts mehr holen, da über 80 Festplatten allein zur Speicherung benötigt werden...
Der Hauptgrund ist jedoch die schlechte Effektivität der Software! Die superschnellen Rechner von heute laufen oft langsamer als die 386er von vor über 20 Jahren!
Die langsamen Programmiersprachen wie im Vergleich
Geschwindigkeitsvergleich Tabelle unten
habe ich hier im Beitrag ja bereits übersprungen. (die würden Wochen brauchen)
Und wenn dann noch
- Fehler in anderen Bibliotheken (gefunden & beseitigt)
- und eine doppelt so schnelle Berechnung wie Mathematica
herausspringen,
hat sich der Aufwand gelohnt.
|
Profil
|
FriedrichLaher
Ehemals Aktiv
Dabei seit: 30.10.2001
Mitteilungen: 1923
Wohnort: Wien,Oesterr., Wohnort Stuttgart
 | Beitrag No.33, eingetragen 2017-12-07
|
@hyperG, #32
\quoteon
Auch habe ich bis heute mein 8stelliges Geburtsdatum darin noch nicht gefunden...
\quoteoff
8 mal je 4Bit
oder 12(Jahr), 4(Monat), 5(Tag)
?
kürzere Folgen haben vermutlich höhere Warscheinlichkeit .
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.34, vom Themenstarter, eingetragen 2017-12-08
|
Datum ist bei mir immer so formatiert:
dd.mm.yyyy also wie auf dieser Internetseite oben bei "Die Mathe-Redaktion" 08.12.2017
Bei Verschlüsselungen des Datums kann man natürlich beliebig weit runtergehen - alles nur eine Frage des Entschlüsselns...
in der anderen Konstante A173201
ist mein Datum bei gleicher Stellenanzahl bereits 4 mal enthalten.
Komisch ist auch, dass bei A133731 bereits ab 100000 Ziffern immer die Ziffer 6
am seltensten vorkommt. Normalerweise ändert sich das ständig...
Ich muss mal eine Analyse der "garantiert bis zur Pos n vorkommenden Stellen" machen, wie das schon bei einigen irrationalen Konstanten zu finden ist:
BisZuWelcherNKalleStringKombi...
Bei Pi braucht man 1816743912 Nachkommastellen, um garantiert ALLE 8stelligen Ziffernkombinationen zu finden.
|
Profil
|
Carmageddon
Senior
Dabei seit: 22.12.2009
Mitteilungen: 664
| Beitrag No.35, eingetragen 2017-12-08
|
\quoteon(2017-12-07 22:22 - hyperG in Beitrag No. 32)
\quoteon(2017-12-07 13:34 - Carmageddon in Beitrag No. 30)
... Warum und für was benötigt man den bitte 200 Millionen Stellen Genauigkeit? ... (reiner Selbstzweck) fällt mir da leider absolut nichts ein? :-?
\quoteoff
Normalerweise rechtfertige ich mich nicht für jeden meiner Rekorde:
133731Ziegenfaktor
Ich muss da mal eine Seite aller Begründungen anlegen, damit ich nicht in jedem Forum oder anderen Treffen immer wieder Fakten neu sammeln muss...
\quoteoff
Niemand verlangt dass du dich rechtfertigst? Ich war einfach nur neugierig - vielleicht steht ja eine konkrete Anwendung dahinter die eben diese Genauigkeit verlangt.
Eine rein intrinsische Motivation ist für mich Antwort genug!
\quoteon(2017-12-07 22:22 - hyperG in Beitrag No. 32)
Und wenn dann noch
- Fehler in anderen Bibliotheken (gefunden & beseitigt)
- und eine doppelt so schnelle Berechnung wie Mathematica
herausspringen,
hat sich der Aufwand gelohnt.
\quoteoff
Das sehe ich ein.
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.36, vom Themenstarter, eingetragen 2017-12-09
|
So, geschafft: mein Geburtstag endlich in den Nachkommastellen des
Ziegenfaktors gefunden!
240 Mio. Dezimalstellen wobei die eine sincos-Berechnung etwa 46 min brauchte.
Ab 250 Mio. die beschriebenen Überläufe - leider immer erst nach 45 min. ganz tief versteckt... das bekomme ich so schnell nicht behoben.
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.37, vom Themenstarter, eingetragen 2017-12-13
|
Da die sin-Obergrenze bei 240 Mio. lag -> weiter mit anderen nichtlinearen Funktionen.
Neue Obergrenze nun 323 Mio. Dezimalstellen: ( Ziegenfaktor )
- [precision is 1072982817 bits -> intern also über 1 Mrd. ! ]
- Validierung mit anderer Iteration & anderen Argumenten
- durch Wegfall der letzten Wandlung {cos(x/2)*2 } weniger als 1h 30 min Gesamtzeit
- Ziffer 8 ist nicht mehr am häufigsten (das war für irrationale Zahlen ja schon nicht mehr normal)
- die zu seltene Ziffer 6 verringerte sich von 0,07406 % Abweichung auf schon fast normale 0,04081%
Interessant: der RAM-Bedarf wächst bei so vielen Nachkommastellen auf fast 16 GB an.
Mit einfachen 8 GB PC nicht mehr machbar!
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.38, vom Themenstarter, eingetragen 2019-01-20
|
Mit mathematica waren 400 Mio. Stellen kein Problem
133731Ziegenfaktor
Größtes Problem war die nachträgliche Vergrößerung der
Nachkommastellen-Anzahl, da der Import immer die Berechnungsgenauigkeit begrenzte.
Vergleiche hier
|
Profil
|
hyperG
Senior
Dabei seit: 03.02.2017
Mitteilungen: 2159
| Beitrag No.39, vom Themenstarter, eingetragen 2019-01-21
|
600 Mio -> erledigt (allerdings nun doch mit mathematica)
1 Mrd -> da reichen meine 32 GB nicht :-(
|
Profil
|
Home / Seite ohne Frame
Wie unterstütze ich den Matheplaneten mit Geld?
