Eigene Funktion mit Exp und log zum Wurzelziehen!

[rustyoldguy]

Hallo Leute!

Kürzlich habe ich mal aus Jux meinen alten Turbo Pascal 7.0 Compiler auf DOSEMU angeworfen.
Dabei stellte ich fest, das dieser keine Eigene Funktion zum Potenzieren wie bei c (pow()) hat.
Wenn man nun etwas herum sieht, so kann man sich selber eine Funktion zum Wurzelziehen mit log()
und exp() basteln.

Falls aber einer mal im Unterricht von seinem Dozenten dann kritisiert wird, das exp und log ja ebenfalls
aus math.h stammen, so sollte in dem Fall dann die Funktionen log() und exp ebenfalls selbst geschrieben
sein:

#include <stdio.h>
#include <cstdlib>
#include <cmath>

double getfakul(double gz);
double radtograd(double boma);
double exponent(double x);
double Power(int n, double x);
int xpow(int x, int y);
double natlogb(double x);
double own_arctanh(double z);
double natlog(double x);
double xwurzely(double x, double y);
void showerge(double x);

double xrooty(double x, double y)
{
double q, a, b;
int c = 0;

c = (int)y % 2;

if(y < 0) b = y * -1;
 else b = y;


q = (1/x) * natlogb(b);
a = exponent(q);

if(c) a *= -1;

return a;
}

//IN(xy) = IN(x) + In(y)
int main(int argc, char **argv)
{
double y,x; //ata, q;

    printf("\nexptest V0.05\n\n");
 
    printf("xpow 10 ^ 4 = %d\n",  xpow(10, 4));
    printf("Power10 ^ 4 = %f\n",  Power(4, 10));
    for (y = 0; y <= 5; y++ )
     printf("-0.6^%f = %.10f\n", y, Power(y, -0.6));
 
  printf("exp5=%f\n", exp(5));
  printf("exp5=%f\n\n", exponent(5));
  x = 0.25;
  showerge(x);
  x = 125;
  showerge(x);
  x = 12;
  showerge(x);
  x = 2;
  showerge(x);
  x = 0.8;
  showerge(x);
  printf("--------------------\n");
  //printf("%f=%f\n", (a-1), natlog(a - 1));
  /*
  q = (0.25-1) / (0.25+1);
  printf("x-1/ x+1=%f\n", q);
  ata = own_arctanh(q);
  printf("In mit own_arctanh 0.25=%f\n", 2*ata);
  */
  printf("3te Wurzel aus 125 mit natlogb und exponent errechnet = %f\n", xwurzely(3, 125));
  printf("3te Wurzel aus 125 mit natlogb und exponent errechnet = %f\n", xrooty(3, -125));
  printf("4te Wurzel aus 125 mit natlogb und exponent errechnet = %f\n", xrooty(4, -125));

 return EXIT_SUCCESS;
}

double getfakul(double gz)
{
double slei, dmg;
dmg = 1;
 for (slei = 1; slei <= gz; slei++)
  dmg *= slei;
return(dmg);
}

// Wandelt Bogenmass in Grad um
double radtograd(double boma)
{
 double grad, pi =M_PI; //= 3.141592653589793;

 grad = ((boma * 180) / pi);
 return grad;
}

double exponent(double x)
{
 double rewer = 0, p,f,i, q;
 
 for (i=0; i <= 50; i++)
  {
    p = Power(i, x);
    f = getfakul(i);
    //printf("%f^%f=%f\n", x,i,p);
    //printf("%f!=%f\n", i,f);
    q = p / f;
    //printf("%f/%f=%f   exp=%f\n", p,f, q, rewer);
 
    rewer += q;
  }
 return rewer;
}
//The exponential value of 5.000000 is 148.413159.

double Power(int n, double x)
{
 const double k=x;
 int MultCount = 0;

 if (n >= 1)
  while (MultCount < n-1)
   {
    MultCount++;
    x *= k;
   }
   else if(n == 0) x = 1;
   else if(n < 0)
    {
     while (MultCount > n + 1)
      {
       MultCount--;
       x *= k;
      //printf("**n=%d    x=%f **\n",n, x);
      }
     x= -1 / x;
    }

    return x;
}

int xpow(int x, int y)
{
 int i, rewer = 1;
 
 
 if (y == 0) return 1;
 if (y == 1) return x;
 
 for (i = 0; i < y; i++)
 {
  rewer *= x;
  printf("x=%d y=%d i=%d  rewer=%d\n", x, y, i,  rewer);
 
 }

return rewer;
}

double natlogb(double x)
{
 double rewer, k, xm1, xp1, p, q, kx;
Diese Funktion wird dann später noch gebraucht, s
 xm1 = x - 1;
 xp1 = x + 1;
 q = xm1 / xp1;
 rewer = q;
 //printf("\nSchleife anfang x=%f\n", x);
 //printf("k=%f q=%f  rewer=%f\n", k, q,rewer);
 for (k = 3; k< 1200;k+= 2)
 {
  p = Power(k, q);
  kx = p * (1/k);
  rewer += kx;
  //printf("k=%f q=%f  p=%.12f  p/k=%.12f rewer=%f\n", k, q,p, kx, rewer);
  }
 //printf("k=%f q=%f  p=%.12f  p/k=%.12f rewer=%f\n", k, q,p, kx, rewer);
 //printf("Schleife ende\n");
 return 2* rewer;
}

double own_arctanh(double z)
{
 double a, i, j, sum = z;  
 
 for (i = 3; i <= 32; i+=2)
  {
   j = Power(i, z);//pow(z, i);
   a = j / i;
   sum += a;
  }

  return sum;
}

double natlog(double x)
{
 double rewer = x, a;
 a = (x-1) / (x+1);
 rewer = 2* own_arctanh(a);

 return rewer;
}

double xwurzely(double x, double y)
{
double q, a;
q = (1/x) * natlogb(y);
a = exponent(q);
return a;
}

void showerge(double x)
 {
  printf("\n\nErgebnis von x:\n");
  printf("log     %f =%.12f\n", x, log(x));
  printf("natlog  %f =%.12f\n", x, natlog(x));
  printf("natlogb %f =%.12f\n", x, natlogb(x));
 }

Nun zur Erklärung der oberen Funktionen:
mit

int xpow(int x, int y)
{
 int i, rewer = 1;

 if (y == 0) return 1;
 if (y == 1) return x;

 for (i = 0; i < y; i++)
 {
  rewer *= x;
 }

return rewer;
}

wäre der einfachste Weg um Zahlen potenzieren zu können.
will man aber damit Wurzel ziehen (2te Wurzel aus x = x hoch (1/2)) ist aber hier schone Essig.

Das geht mit exp und log besser:

double xwurzely(double x, double y)
{
double q, a;
q = (1/x) * log(y);
a = exp(q);
return a;
}

Wohl aber kann hier der Dozent bemängeln, das exp und log ebenfalls nicht selbst geschrieben seien.
Hier also zunächst die Funktion exp(onent) in eigener Fassung:

double exponent(double x)
{
 double rewer = 0, p,f,i, q;
 
 for (i=0; i <= 50; i++)
  {
    p = Power(i, x);
    f = getfakul(i);
    //printf("%f^%f=%f\n", x,i,p);
    //printf("%f!=%f\n", i,f);
    q = p / f;
    //printf("%f/%f=%f   exp=%f\n", p,f, q, rewer);
 
    rewer += q;
  }
 return rewer;
}

Aber es wird noch eine eigene Fassung von log, der Funktion zum Berechnen des natürlichen Logarithmus benötigt.
Dazu sind im Beispiel zwei Fassungen vorhanden, natlog und natlogb.
Zuerst die Fassung mit natlog, welche den arctanh benötigt:

double own_arctanh(double z)
{
 double a, i, j, sum = z;  
 
 for (i = 3; i <= 32; i+=2)
  {
   j = Power(i, z);//pow(z, i);
   a = j / i;
   sum += a;
  }

  return sum;
}

double natlog(double x)
{
 double rewer = x, a;
 a = (x-1) / (x+1);
 rewer = 2* own_arctanh(a);

 return rewer;
}

Dann noch die Funktion natlogb:

double natlogb(double x)
{
 double rewer, k, xm1, xp1, p, q, kx;
 
 xm1 = x - 1;
 xp1 = x + 1;
 q = xm1 / xp1;
 rewer = q;
 //printf("\nSchleife anfang x=%f\n", x);
 //printf("k=%f q=%f  rewer=%f\n", k, q,rewer);
 for (k = 3; k< 1200;k+= 2)
 {
  p = Power(k, q);
  kx = p * (1/k);
  rewer += kx;
  //printf("k=%f q=%f  p=%.12f  p/k=%.12f rewer=%f\n", k, q,p, kx, rewer);
  }
 //printf("k=%f q=%f  p=%.12f  p/k=%.12f rewer=%f\n", k, q,p, kx, rewer);
 //printf("Schleife ende\n");
 return 2* rewer;
}

Bei der Funktion Power ist zu beachten, das bei meiner Version die Parameter bei der Übergabe vertauscht sind.
x hoch y wäre hier Power(y,x);

Die zweite Version einer Funktion zum Wurzelziehen berücksichtigt, das bei einer negativen Zahl
die Wurzel, falls der Exponent ungerade ist, der Wurzelwert ebenfalls negativ ist.
Beispiel:

(-5)^2 = (+25)
(-5)^3 = (-125)

Die dritte Wurzel aus (-125) erhält man mit (-125)^(1/3). Das ist dann (-5).

Folgende Version verbessert hier die Ausgabe:

double xrooty(double x, double y)
{
double q, a, b;
int c = 0;

c = (int)y % 2;

printf("c: %d\n", c);

if(y < 0) b = y * -1;
 else b = y;


q = (1/x) * natlogb(b);
a = exponent(q);

if(c < 0) a *= -1;

return a;
}
}

Bei mir am Terminal sieht das Ergebnis dann so aus:
exptest V0.05

x=10 y=4 i=0 rewer=10
x=10 y=4 i=1 rewer=100
x=10 y=4 i=2 rewer=1000
x=10 y=4 i=3 rewer=10000
xpow 10 ^ 4 = 10000
Power10 ^ 4 = 10000.000000
-0.6^0.000000 = 1.0000000000
-0.6^1.000000 = -0.6000000000
-0.6^2.000000 = 0.3600000000
-0.6^3.000000 = -0.2160000000
-0.6^4.000000 = 0.1296000000
-0.6^5.000000 = -0.0777600000
exp5=148.413159
exp5=148.413159



Ergebnis von x:
log 0.250000 =-1.386294361120
natlog 0.250000 =-1.386294356735
natlogb 0.250000 =-1.386294361120


Ergebnis von x:
log 125.000000 =4.828313737302
natlog 125.000000 =4.282999474420
natlogb 125.000000 =4.828313737075


Ergebnis von x:
log 12.000000 =2.484906649788
natlog 12.000000 =2.484144862627
natlogb 12.000000 =2.484906649788


Ergebnis von x:
log 2.000000 =0.693147180560
natlog 2.000000 =0.693147180560
natlogb 2.000000 =0.693147180560


Ergebnis von x:
log 0.800000 =-0.223143551314
natlog 0.800000 =-0.223143551314
natlogb 0.800000 =-0.223143551314
--------------------
3te Wurzel aus 125 mit natlogb und exponent errechnet = 5.000000
3te Wurzel aus 125 mit natlogb und exponent errechnet = -5.000000
4te Wurzel aus 125 mit natlogb und exponent errechnet = -3.343702

Damit dürte dann Kritikern, welche bemängeln, das einige Funktionen ebenfalls aus math.h stammen,
der Wind aus den Segeln genommen sein. Zumindest sieht dann diese Version eines Source-Codes besser
aus als wie mit der einfachen Funktion xpow. Die zweite Version hat aber einen Fehler:
Bei einem geradzahligen Exponenten wären bei positiven Potenzen negative sowohl positive Wurzelwerte möglich.
Aber die Standardfunktion sqrt() von Compilern berücksichtigen dies ebenfalls nicht. Man sollte darüber Bescheid wissen.

Würde man die erste Version mit der letzten tauschen und diese Zeilen im Beispiel implementieren, so erhielte man:
3te Wurzel aus +125 mit natlogb und exponent errechnet = 5.000000
c: -1
3te Wurzel aus -125 mit natlogb und exponent errechnet = -5.000000
c: 1
4te Wurzel aus +125 mit natlogb und exponent errechnet = 3.343702
c: 1
5te Wurzel aus +625 mit natlogb und exponent errechnet = 3.620573
c: -1
5te Wurzel aus -625 mit natlogb und exponent errechnet = -3.620573
Hit any key to continue...


Übrigens:

Ich hatte die Disketten in Keksdosen aufbewahrt. Diese sind aus Metall.
Darin befanden sich auch die Startdisketten zu Olivetti-Unix 2.4.1 von 1991.
Mit Hilfe von Linux habe ich daraus je eine Image-Datei gemacht.
Das Band zum Betriebssystem habe ich leider nicht mehr, ging verloren.
Aber zum Testen langt es aus. Leider sind diese vier Disketten selbst gepackt als Anhang zu groß. Warum auch immer.
Hier die zwei Startdisketten für die AT-Version. Die SCSI-Version habe ich ebenfalls zu einer Image-Datei kopiert.
Stelle ich aber nur auf Anfrage hier rein.

Ich möchte keinesfalls Dozenten in die Suppe spucken. Nur dauert es eben, bis man die
richtigen Lösungen im Web gefunden hat. Das geht von der Zeit zum Programmieren weg.
Meine Funktionen haben auch Fehler. Aber wer diese in seine Sourcen einbaut, der wird diese
sowieso extra testen.

Das Ganze stelle ich unter GPL-Lizenz


Viel Spaß wünscht euch
rustyoldguy