diff src/decimal/math.d @ 4:b37c218c1442

Initial development of BCD integers
author Paul (paul.d.anderson@comcast.net)
date Thu, 18 Mar 2010 18:10:25 -0700
parents a984d3056cc4
children
line wrap: on
line diff
--- a/src/decimal/math.d	Thu Mar 18 18:09:20 2010 -0700
+++ b/src/decimal/math.d	Thu Mar 18 18:10:25 2010 -0700
@@ -34,6 +34,7 @@
 
 import decimal.decimal;
 import decimal.context;
+//import std.math;
 
 public:
     
@@ -42,6 +43,14 @@
     // CONSTANTS
     //
     //--------------------------------
+	
+/+	public Decimal HALF;
+	public Decimal ONE = Decimal(1);
+	public Decimal ZERO = Decimal(0);
+	
+	static {
+		HALF = Decimal("0.5");
+	}+/
 
     /**
      * Returns the value of e to the default precision.
@@ -59,14 +68,62 @@
         return result;
     }
      
+/+
     /**
      * Returns the value of pi to the default precision.
      */    
     Decimal pi() {
-        Decimal result;
-        return result;
+		Decimal ONE = Decimal(1);
+		writeln("1 = ", ONE);
+		Decimal TWO = Decimal(2);
+		writeln("2 = ", TWO);
+		Decimal HALF = Decimal(0.5);
+		writeln("0.5 = ", HALF);
+        Decimal x = sqrt(TWO);
+		writeln("x = ", x);
+		Decimal y = sqrt(x);
+		writeln("y = ", y);
+		Decimal p = TWO + x;
+		writeln("p = ", p);
+		x = y;
+		int i = 0;
+		while (true) {
+			writeln("i = ", i);
+			x = HALF * (x + ONE/x);
+			writeln("step 1");
+			Decimal np = p * ((x + ONE)/(y + ONE));
+			writeln("step 2");
+			writeln("np = ", np);
+			if (p == np) return p;
+			writeln("step 3");
+			p = np;
+			writeln("step 4");
+			x = sqrt(x);
+			writeln("step 5");
+//			writeln("x = ", x);
+			Decimal oox = ONE/x;
+//			writeln("ONE/x = ", oox);
+//			writeln ("x + ONE/x = ", x + oox);
+//			Decimal t1 = oox + x;
+//			writeln("t1 = ", t1);
+//			Decimal t1 = x + oox; // ONE + oox; //x + x; // + ONE/x;
+			writeln("step 6");
+//			Decimal t2 = (y  * x) + ONE; ///x; //ONE / (y + ONE);
+//			y = ONE/x + y * x;
+			writeln("step 7");
+			y = (ONE/x + y * x) / (y + ONE); //t1 / t2;
+			writeln("step 8");
+			i++;
+//			break;
+		}
+        return p;
     }
-     
+	
+	unittest {
+		write("pi....");
+		writeln("pi = ", pi());
+	}
++/     
     /**
      * Returns the value of pi to the specified precision.
      */    
@@ -76,19 +133,38 @@
     }
      
     /**
-     * Returns the square root of the argument to the default precision.
+     * Returns the square root of the argument to the current precision.
      */    
     Decimal sqrt(Decimal arg) {
-        Decimal result;
-        return result;
+		return sqrt(arg, context.precision);
     }
+	
+	
+	unittest {
+		write("sqrt.....");
+		Decimal dcm = Decimal(4);
+		assert(sqrt(dcm) == Decimal(2));
+		writeln("sqrt(2) = ", sqrt(Decimal(2)));
+	}
      
     /**
      * Returns the square root of the argument to the specified precision.
      */    
     Decimal sqrt(Decimal arg, uint precision) {
-        Decimal result;
-        return result;
+		pushContext();
+		context.precision = precision + 2;
+		Decimal HALF = Decimal(0.5);
+		Decimal ONE = Decimal(1);
+		Decimal x = HALF * (arg + ONE);
+		Decimal xp;
+		while(true){
+			xp = x;
+			x = HALF * (x + arg/x);
+			if (x == xp) break;
+		}
+		popContext();
+		round(xp);
+        return xp;
     }
     
     //--------------------------------
@@ -312,7 +388,7 @@
         return result;
     }
 
-    //--------------------------------
+/+    //--------------------------------
     //
     // General Decimal Arithmetic Specification Functions
     //
@@ -553,6 +629,6 @@
         Decimal result;
         return result;
     }
++/
 
 
-