From 04e68443040c7abad84d66477e98f93bed701760 Mon Sep 17 00:00:00 2001 From: Lars-Dominik Braun Date: Mon, 4 Feb 2019 13:09:03 +0100 Subject: Initial import --- system/base/1.7.5/src/mathlib | 268 ++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 268 insertions(+) create mode 100644 system/base/1.7.5/src/mathlib (limited to 'system/base/1.7.5/src/mathlib') diff --git a/system/base/1.7.5/src/mathlib b/system/base/1.7.5/src/mathlib new file mode 100644 index 0000000..c726495 --- /dev/null +++ b/system/base/1.7.5/src/mathlib @@ -0,0 +1,268 @@ +(* ------------------- VERSION 2 06.03.86 ------------------- *) +PACKET mathlib DEFINES sqrt, **, exp, ln, log2, log10, e, pi, + sin, cos, tan, sind, cosd, tand, + arctan, arctand, random, initializerandom : + +LET pii = 3.141592653589793238462, + pi2 = 1.570796326794896619231, + pi3 = 1.047197551196597746154, + pi6 = 0.523598775598298873077, + pi4 = 1.273239544735162686151, + ln2 = 0.693147180559945309417, + lg2 = 0.301029995663981195213, + ln10 = 2.302585092994045684018, + lge = 0.434294481903251827651, + ei = 2.718281828459045235360, + pi180 = 57.295779513082320876798, + sqrt3 = 1.732050807568877293527, + sqr3 = 0.577350269189625764509, + sqr3p2= 3.732050807568877293527, + sqr3m2= 0.267949192431122706473, + sqr2 = 0.707106781186547524400; + +REAL VAR rdg::0.4711; + +REAL PROC pi: pii END PROC pi; +REAL PROC e : ei END PROC e; + +REAL PROC ln ( REAL CONST x ): + log2(x) * ln2 +END PROC ln; + +REAL PROC log10( REAL CONST x ): + log2(x) * lg2 +END PROC log10; + +REAL PROC log2 ( REAL CONST z ): + REAL VAR t, summe::0.0, x::z; + IF x=1.0 THEN 0.0 + ELIF x>0.0 THEN normal + ELSE errorstop("log2: " + text (x,20)); 0.0 FI. + +normal: + IF x >= 0.5 THEN normalise downwards + ELSE normalise upwards FI; + IF x < sqr2 THEN summe := summe - 0.75; t := trans8 + ELSE summe := summe - 0.25; t := trans2 FI; + summe + reihenentwicklung. + + normalise downwards: + WHILE x >= 8.0 REP x := 0.0625 * x; summe:=summe+4.0 PER; + WHILE x >= 1.0 REP x := 0.5 * x; summe:=summe+1.0 PER. + + normalise upwards: + WHILE x<=0.0625 REP x := 16.0 * x; summe:=summe-4.0 PER; + WHILE x<= 0.5 REP x := 2.0 * x; summe:=summe-1.0 PER. + + trans8: (x - 0.5946035575013605)/(x + 0.5946035575013605). + trans2: (x - 0.8408964152537145)/(x + 0.8408964152537145). + + reihenentwicklung: x := t * t; t * 0.06405572387119384648 * + ((((((3.465*x+4.095)*x+5.005)*x+6.435)*x+9.009)*x+15.015)*x+45.045) +END PROC log2; + +REAL PROC sqrt ( REAL CONST z ): + REAL VAR y0, y1, x::z; + INT VAR p :: decimal exponent(x) DIV 2; + IF p <= -64 THEN 0.0 + ELIF x < 0.0 THEN errorstop("sqrt: " + text (x,20)); 0.0 + ELSE nontrivial FI. + + nontrivial: + set exp (decimal exponent (x) -p-p, x); + IF x<10.0 THEN x := 5.3176703 - 40.760905/( 8.408065 + x ) + ELSE x := 16.81595 - 1288.973 /( 84.08065 + x ) FI; + y0 := x; + set exp (decimal exponent (x) + p, y0); + y1 := 0.5 * ( y0 + z/y0 ); + y0 := 0.5 * ( y1 + z/y1 ); + y1 := 0.5 * ( y0 + z/y0 ); + 0.5 * ( y1 + z/y1 ) +END PROC sqrt; + +REAL PROC exp ( REAL CONST z ): + REAL VAR x::z, a::1.0; BOOL VAR negativ :: x<0.0; + IF negativ THEN x := -x FI; + IF x>292.42830676 + THEN IF NOT negativ THEN errorstop ("REAL-Ueberlauf") FI ; 0.0 + ELIF x<=0.0001 + THEN ( 0.5*z + 1.0 ) * z + 1.0 + ELSE approx + FI. + + approx: + IF x > ln10 + THEN x := lge*x; + a := 1.0; + set exp (int(x), a); + x := frac(x)*ln10 + FI; + IF x >= 2.0 THEN a := 7.389056098930650227230*a; x := x-2.0 FI; + IF x >= 1.0 THEN a := 2.718281828459045235360*a; x := x-1.0 FI; + IF x >= 0.5 THEN a := 1.648721270700128146848*a; x := x-0.5 FI; + IF x >= 0.25 THEN a := 1.284025416687741484073*a; x := x-0.25 FI; + IF x >= 0.125 THEN a := 1.133148453066826316829*a; x := x-0.125 FI; + IF x >= 0.0625THEN a := 1.064494458917859429563*a; x := x-0.0625FI; + a:=a/50.4*(((((((0.01*x+0.07)*x+0.42)*x+2.1)*x+8.4)*x+25.2)*x+50.4)*x+50.4); + IF negativ THEN 1.0/a ELSE a FI . + +ENDPROC exp ; + +REAL PROC tan (REAL CONST x): + IF x < 0.0 THEN - tg( -x * pi4) + ELSE tg( x * pi4) FI +END PROC tan; + +REAL PROC tand (REAL CONST x): + IF x < 0.0 THEN - tg( -x / 45.0) + ELSE tg( x / 45.0) FI +END PROC tand; + +REAL PROC tg (REAL CONST x ): + REAL VAR q::floor(x), s::x-q; INT VAR n; + q := q - floor(0.25*q) * 4.0 ; + IF q < 2.0 + THEN IF q < 1.0 + THEN n:=0; + ELSE n:=1; s := 1.0 - s FI + ELSE IF q < 3.0 + THEN n:=2; + ELSE n:=3; s := 1.0 - s FI + FI; + q := s * s; + q := (((((((((-5.116186989653120e-11*q-5.608325022830701e-10)*q- + 9.526170109403018e-9)*q-1.517906721393745e-7)*q-2.430939946375515e-6)*q- + 3.901461426385464e-5)*q-6.324811612385572e-4)*q-1.076606829172646e-2)*q- + 0.2617993877991508)*q+pi4); + + SELECT n OF + CASE 0 : s/q + CASE 1 : q/s + CASE 2 : -q/s + OTHERWISE : -s/q ENDSELECT . + +END PROC tg; + +REAL PROC sin ( REAL CONST x ): + REAL VAR y, r, q; + IF x < 0.0 THEN y := -x; q := 4.0 ELSE y := x; q := 0.0 FI; + y := y * pi4; + r := floor(y); + sincos( q+r , y-r ) +END PROC sin; + +REAL PROC sind ( REAL CONST x ): + REAL VAR y, r, q; + IF x < 0.0 THEN y := -x; q := 4.0 ELSE y := x; q := 0.0 FI; + y := y / 45.0; + r := floor(y); + sincos( q+r , y-r ) +END PROC sind; + +REAL PROC cos ( REAL CONST x ): + REAL VAR y, q; + IF x < 0.0 THEN y := -x ELSE y := x FI; + y := y * pi4; + q := floor(y); + sincos( q+2.0, y-q ) +END PROC cos; + +REAL PROC cosd ( REAL CONST x ): + REAL VAR y, q; + IF x < 0.0 THEN y := -x ELSE y := x FI; + y := y / 45.0; + q := floor(y); + sincos( q+2.0, y-q ) +END PROC cosd; + +REAL PROC sincos ( REAL CONST q, y ): + REAL VAR r :: q - floor( 0.125*q + 0.1 ) * 8.0; + IF r >= 4.0 THEN IF r >= 6.0 THEN IF r >= 7.0 THEN - sin approx(1.0-y) + ELSE - cos approx(y) FI + ELSE IF r >= 5.0 THEN - cos approx(1.0-y) + ELSE - sin approx(y) FI FI + ELSE IF r >= 2.0 THEN IF r >= 3.0 THEN sin approx(1.0-y) + ELSE cos approx(y) FI + ELSE IF r >= 1.0 THEN cos approx(1.0-y) + ELSE sin approx(y) FI FI FI +END PROC sincos; + +REAL PROC sin approx ( REAL CONST x ): + REAL VAR z::x*x; + x*((((((0.6877101540593035e-11*z-0.1757149296873372e-8)*z+0.3133616216672568 + e-6)*z-0.3657620415845891e-4)*z+0.2490394570188737e-2)*z-0.807455121882e-1)* + z+0.7853981633974483) +END PROC sin approx; + +REAL PROC cos approx ( REAL CONST x ): + REAL VAR z::x*x; + ((((((-0.3857761864560276e-12*z+0.115004970178141e-9)*z-0.246113638267419e-7 + )*z+0.3590860445885748e-5)*z-0.3259918869266875e-3)*z+0.1585434424381541e-1) + *z-0.3084251375340425)*z+1.0 +END PROC cos approx; + +REAL PROC arctan ( REAL CONST y ): + REAL VAR f, z, x; BOOL VAR neg :: y < 0.0; + IF neg THEN x := -y ELSE x := y FI; + IF x>1.0 THEN f := a ELSE f := -b; neg := NOT neg FI; + z := x * x; + x := x/(((((((0.0107090276046822*z-0.01647757182108040)*z + +0.02177846332482151)*z-0.03019339673273880)*z+0.04656083561183398)*z + -0.0888888888888888)*z+0.3333333333333333)*z+1.0); + IF neg THEN x - f ELSE f - x FI. + + a:IF x>sqr3p2 THEN x := 1.0/x; pi2 ELSE x := 4.0/(sqrt3+x+x+x)-sqr3; pi3 FI. + b:IF x0 REP + m := n DIV 2 ; + IF m + m = n + THEN n := m ; + r := r*r + ELSE n DECR 1 ; + p := p*r + FI + END REP ; + IF b>0 + THEN p + ELSE 1.0 / p + FI + FI . + +END OP ** ; + +REAL PROC random: + rdg:=rdg+pii;rdg:=rdg*rdg;rdg:=rdg*rdg;rdg:=rdg*rdg;rdg:=frac(rdg);rdg +END PROC random; + +PROC initializerandom ( REAL CONST z ): + rdg := frac(z) +END PROC initializerandom; + +END PACKET mathlib; + -- cgit v1.2.3