sqrt funktion nachbauen



  • hi, welche methode ist schneller um die wurzel zu berechnen?


  • Mod

    algoman1 schrieb:

    hi, welche methode ist schneller um die wurzel zu berechnen?

    Erst einmal Newtons Methode probieren, damit liegt man sicherlich nicht schlecht. Das eingebaute sqrt dürfte aber schneller sein. Die benutzen zwar wahrscheinlich auch eine Variation von Newtons Methode, aber noch zusätzlich irgendwelche Mikrooptimierungen. Oder direkt Hardwareroutinen der CPU, die zwar auch ähnlich funktionieren sollten, aber natürlich nochmal viel schneller sind.



  • Es wird häufig auch CORDIC eingesetzt.



  • hi, eine moeglichkeit ist die newton raphson's method...
    warum wird am anfang x0 = a/2 zugewiesen?

    double sqrt(double a) {
       double x0 = a/2;
       double err = 1, eps = Math.pow(10,-11);
       while (err >= eps) {
           double x = x0 - ((x0 * x0 - a)/(2 * x0));
           err = x-x0;
           x0 = x;
       }
       return x0;
    }
    


  • Ehrlich gesagt ist das Scheißegal, weil das Teil sowieso ab ca.der dritten Nachkommastelle falsche Ergebnisse liefert:

    2: 1.41666666666666674068
    

    mit der obigen Funktion

    2: 1.41421356237309514547
    

    mit der C-Funktion sqrt

    5: 2.25000000000000000000
    

    mit der obigen Funktion

    5: 2.23606797749978980505
    

    mit der C-Funktion sqrt

    Das scheinen alles, inklusive Newtons Formel, mehr oder weniger Variationen der Iterationsformel von Archimedes zu sein. Probiers doch einfach mal mit dem Original 😉

    Viele Grüße

    knoppi



  • knoppi schrieb:

    Ehrlich gesagt ist das Scheißegal, weil das Teil sowieso ab ca.der dritten Nachkommastelle falsche Ergebnisse liefert:

    Ob das an der 'percision' liegen könnte? :p

    Das scheinen alles, inklusive Newtons Formel, mehr oder weniger Variationen der Iterationsformel von Archimedes zu sein.

    Die Archimedes-Formel kenn ich zwar nicht, aber vielleicht meinst du das Heron-Verfahren. Das ist in der Tat ein Spezialfall des Newton-Verfahrens.



  • Iterationsformel von Archimedes von Syrakus

    xneu = 0.5*(xalt+a/xalt)
    
    #include <math.h>
    #include <stdio.h>
    
    /* Errechnen der Quadradwurzel nach der Iterationsformel von Archimedes */
    
    double archimed(double);
    
    int main(int argc, char *argv[]) {
    	double a = 0;
    
    	printf("Von welcher Zahl soll die Quadratwurzel gezogen werden ?\n");
    	scanf("%lf", &a);
    
    	printf("\nDie Quadratwurzel von %.2lf ist %.11lf\n", a, archimed(a));
    	printf("\nDie Quadratwurzel von %.2lf ist %.11lf\n", a, sqrt(a));
    
    	return 0;
    }
    
    double archimed(double a) {
    	double xneu = 0, xalt = 0, epsilon = 0.00000000001;
    
    	xneu = a;
    	do {
    		xalt = xneu;
    		xneu = 0.5*(xalt+a/xalt);
    	} while((xalt-xneu) > epsilon);
    
    	return xneu;
    }
    

    Ich setze hier mal eine Eingabe für a >= 1 voraus, um grad das Wesentliche zu zeigen.

    Ich hab den Vergleich mit sqrt() mit einbezogen. Außerdem habe ich noch Vergleiche mit Taschenrechner und dem hier auf meinem PC gemacht. Ist zwar nicht unbedingt ein Beweis, aber bisher stimmte alles.

    Viele Grüße

    knoppi


  • Mod

    Was ist, wenn ich die Wurzel von 10^-50 bestimmen will?



  • knoppi schrieb:

    Iterationsformel von Archimedes von Syrakus

    xneu = 0.5*(xalt+a/xalt)
    

    Auch bekannt als Heron-Verfahren. Hatte ich glaub ich schon erwähnt... manche nennen es auch babylonisches Verfahren.


  • Mod

    SeppJ schrieb:

    Was ist, wenn ich die Wurzel von 10^-50 bestimmen will?

    Um mal das Problem zu verdeutlichen habe ich ein kleines Testprogramm geschrieben. Wenn wir die Standardwurzel mal als exakt ansehen, dann vergleichen wir mal Ergebnisse. Eine absolute Präzisionsvorgabe ist schlecht.

    Als Vergleich habe ich mal die berühmt berüchtigte Q3-Wurzel genommen (auch ein Newtonverfahren) mit einmal 1 Schritt und einmal 2 Schritten, Schrittzahl fest, ohne feste Präzisionsvorgabe. Man beachte die absoluten und relativen Abweichungen von der Standardwurzel der verschiedenen Verfahren für verschiedene Größenordnungen der Eingabe:

    #include <math.h>
    #include <stdio.h>
    #include <stdlib.h>
    #include <stdint.h>
    
    double sqrt_knoppi(double a) {
        double xneu = 0, xalt = 0, epsilon = 0.00000000001;
    
        xneu = a;
        do {
            xalt = xneu;
            xneu = 0.5*(xalt+a/xalt);
        } while((xalt-xneu) > epsilon);
    
        return xneu;
    }
    
    const int64_t magic_number = 0x5fe6eb50c7b537a9;
    double  sqrt_q3(const double x)
    {
      const double xhalf = 0.5f*x; 
      union
      {
        double x;
        int64_t i;
      } u;
      u.x = x;
      u.i = magic_number - (u.i >> 1); 
      return x*u.x*(1.5f - xhalf*u.x*u.x);
    }   
    
    double  sqrt_q3_2(const double x)
    {
      const double xhalf = 0.5f*x; 
      union
      {
        double x;
        int64_t i;
      } u;
      u.x = x;
      u.i = magic_number - (u.i >> 1); 
      u.x = u.x*(1.5f - xhalf*u.x*u.x); // Mehr von dieser Zeile für höhere Präzision
      return x*u.x*(1.5f - xhalf*u.x*u.x);
    }   
    
    void print_diff(double result, double sqrt)
    {
      printf("%.15f \tdiff = %.15f\tRel = %.15f\n", result, fabs(sqrt -result), fabs(sqrt -result)/sqrt);
    }
    
    void sqrt_test(double d)
    {
      double result = sqrt(d);
      printf("sqrt(%.15f) = %.15f\n", d, result);
      printf("knoppi: sqrt(%.15f) = ", d); print_diff(sqrt_knoppi(d), result);
      printf("knoppi: sqrt_q3(%.15f) = ", d); print_diff(sqrt_q3(d), result);
      printf("knoppi: sqrt_q3_2(%.15f) = ", d); print_diff(sqrt_q3_2(d), result);
      putchar('\n');
    }
    
    int main()
    {
      double d = 15;
      for(int i = 0; i < 30; ++i)
        {
          sqrt_test(d);
          d /= 3;      
        }
    }
    

    https://ideone.com/X65xgk


Anmelden zum Antworten