You are viewing a plain text version of this content. The canonical link for it is here.
Posted to mapreduce-commits@hadoop.apache.org by sz...@apache.org on 2009/07/01 03:35:56 UTC

svn commit: r790015 [2/2] - in /hadoop/mapreduce/trunk: ./ src/examples/org/apache/hadoop/examples/ src/examples/org/apache/hadoop/examples/pi/ src/examples/org/apache/hadoop/examples/pi/math/ src/test/mapred/org/apache/hadoop/examples/pi/ src/test/map...

Added: hadoop/mapreduce/trunk/src/examples/org/apache/hadoop/examples/pi/math/Montgomery.java
URL: http://svn.apache.org/viewvc/hadoop/mapreduce/trunk/src/examples/org/apache/hadoop/examples/pi/math/Montgomery.java?rev=790015&view=auto
==============================================================================
--- hadoop/mapreduce/trunk/src/examples/org/apache/hadoop/examples/pi/math/Montgomery.java (added)
+++ hadoop/mapreduce/trunk/src/examples/org/apache/hadoop/examples/pi/math/Montgomery.java Wed Jul  1 01:35:55 2009
@@ -0,0 +1,78 @@
+/**
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements.  See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership.  The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License.  You may obtain a copy of the License at
+ *
+ *     http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package org.apache.hadoop.examples.pi.math;
+
+/** Montgomery method.
+ * 
+ * References:
+ * 
+ * [1] Richard Crandall and Carl Pomerance.  Prime Numbers: A Computational 
+ *     Perspective.  Springer-Verlag, 2001. 
+ * 
+ * [2] Peter Montgomery.  Modular multiplication without trial division.
+ *     Math. Comp., 44:519-521, 1985.
+ */
+class Montgomery {
+  protected final Product product = new Product();
+
+  protected long N;
+  protected long N_I;  // N'
+  protected long R;
+  protected long R_1;  // R - 1
+  protected int s;
+
+  /** Set the modular and initialize this object. */
+  Montgomery set(long n) {
+    if (n % 2 != 1)
+      throw new IllegalArgumentException("n % 2 != 1, n=" + n);
+    N = n;
+    R = Long.highestOneBit(n) << 1;
+    N_I = R - Modular.modInverse(N, R);
+    R_1 = R - 1;
+    s = Long.numberOfTrailingZeros(R);
+    return this;
+  }
+
+  /** Compute 2^y mod N for N odd. */
+  long mod(final long y) {
+    long p = R - N; 
+    long x = p << 1;
+    if (x >= N) x -= N;
+    
+    for(long mask = Long.highestOneBit(y); mask > 0; mask >>>= 1) {
+      p = product.m(p, p);
+      if ((mask & y) != 0) p = product.m(p, x);
+    }
+    return product.m(p, 1);
+  }
+
+  class Product {
+    private final LongLong x = new LongLong();
+    private final LongLong xN_I = new LongLong();
+    private final LongLong aN = new LongLong();
+
+    long m(final long c, final long d) {
+      LongLong.multiplication(x, c, d);
+      // a = (x * N')&(R - 1) = ((x & R_1) * N') & R_1
+      final long a = LongLong.multiplication(xN_I, x.and(R_1), N_I).and(R_1);
+      LongLong.multiplication(aN, a, N);
+      final long z = aN.plusEqual(x).shiftRight(s);
+      return z < N? z: z - N;      
+    }
+  }
+}
\ No newline at end of file

Added: hadoop/mapreduce/trunk/src/examples/org/apache/hadoop/examples/pi/math/Summation.java
URL: http://svn.apache.org/viewvc/hadoop/mapreduce/trunk/src/examples/org/apache/hadoop/examples/pi/math/Summation.java?rev=790015&view=auto
==============================================================================
--- hadoop/mapreduce/trunk/src/examples/org/apache/hadoop/examples/pi/math/Summation.java (added)
+++ hadoop/mapreduce/trunk/src/examples/org/apache/hadoop/examples/pi/math/Summation.java Wed Jul  1 01:35:55 2009
@@ -0,0 +1,242 @@
+/**
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements.  See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership.  The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License.  You may obtain a copy of the License at
+ *
+ *     http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package org.apache.hadoop.examples.pi.math;
+
+import java.util.ArrayList;
+import java.util.List;
+
+import org.apache.hadoop.examples.pi.Combinable;
+import org.apache.hadoop.examples.pi.Container;
+import org.apache.hadoop.examples.pi.Util;
+
+/** Represent the summation \sum \frac{2^e \mod n}{n}. */
+public class Summation implements Container<Summation>, Combinable<Summation> {
+  /** Variable n in the summation. */
+  public final ArithmeticProgression N;
+  /** Variable e in the summation. */
+  public final ArithmeticProgression E;
+  private Double value = null;
+
+  /** Constructor */
+  public Summation(ArithmeticProgression N, ArithmeticProgression E) {
+    if (N.getSteps() != E.getSteps()) {
+      throw new IllegalArgumentException("N.getSteps() != E.getSteps(),"
+          + "\n  N.getSteps()=" + N.getSteps() + ", N=" + N
+          + "\n  E.getSteps()=" + E.getSteps() + ", E=" + E);
+    }
+    this.N = N;
+    this.E = E;
+  }
+
+  /** Constructor */
+  Summation(long valueN, long deltaN, 
+            long valueE, long deltaE, long limitE) {
+    this(valueN, deltaN, valueN - deltaN*((valueE - limitE)/deltaE),
+         valueE, deltaE, limitE);
+  }
+
+  /** Constructor */
+  Summation(long valueN, long deltaN, long limitN, 
+            long valueE, long deltaE, long limitE) {
+    this(new ArithmeticProgression('n', valueN, deltaN, limitN),
+         new ArithmeticProgression('e', valueE, deltaE, limitE));
+  }
+
+  /** {@inheritDoc} */
+  @Override
+  public Summation getElement() {return this;}
+
+  /** Return the number of steps of this summation */
+  long getSteps() {return E.getSteps();}
+  
+  /** Return the value of this summation */
+  public Double getValue() {return value;}
+  /** Set the value of this summation */
+  public void setValue(double v) {this.value = v;}
+
+  /** {@inheritDoc} */
+  @Override
+  public String toString() {
+    return "[" + N + "; " + E + (value == null? "]": "]value=" + Double.doubleToLongBits(value)); 
+  }
+
+  /** {@inheritDoc} */
+  @Override
+  public boolean equals(Object obj) {
+    if (obj == this)
+      return true;
+    if (obj != null && obj instanceof Summation) {
+      final Summation that = (Summation)obj;
+      return this.N.equals(that.N) && this.E.equals(that.E);
+    }
+    throw new IllegalArgumentException(obj == null? "obj == null":
+      "obj.getClass()=" + obj.getClass());
+  }
+
+  /** Not supported */
+  @Override
+  public int hashCode() {
+    throw new UnsupportedOperationException();
+  }
+
+  /** Covert a String to a Summation. */
+  public static Summation valueOf(final String s) {
+    int i = 1;
+    int j = s.indexOf("; ", i);
+    if (j < 0)
+      throw new IllegalArgumentException("i=" + i + ", j=" + j + " < 0, s=" + s);
+    final ArithmeticProgression N = ArithmeticProgression.valueOf(s.substring(i, j));
+
+    i = j + 2;
+    j = s.indexOf("]", i);
+    if (j < 0)
+      throw new IllegalArgumentException("i=" + i + ", j=" + j + " < 0, s=" + s);
+    final ArithmeticProgression E = ArithmeticProgression.valueOf(s.substring(i, j));
+
+    final Summation sigma = new Summation(N, E);
+    i = j + 1;
+    if (s.length() > i) {
+      final String value = Util.parseStringVariable("value", s.substring(i));
+      sigma.setValue(value.indexOf('.') < 0?
+          Double.longBitsToDouble(Long.parseLong(value)):
+          Double.parseDouble(value));
+    }
+    return sigma;
+  }
+
+  /** Compute the value of the summation. */
+  public double compute() {
+    if (value == null)
+      value = N.limit <= MAX_MODULAR? compute_modular(): compute_montgomery();
+    return value;
+  }
+
+  private static final long MAX_MODULAR = 1L << 32;
+  /** Compute the value using {@link Modular#mod(long, long)}. */
+  double compute_modular() {
+    long e = E.value;
+    long n = N.value;
+    double s = 0;
+    for(; e > E.limit; e += E.delta) {
+      s = Modular.addMod(s, Modular.mod(e, n)/(double)n);
+      n += N.delta;
+    }
+    return s;
+  }
+
+  final Montgomery montgomery = new Montgomery();
+  /** Compute the value using {@link Montgomery#mod(long)}. */
+  double compute_montgomery() {
+    long e = E.value;
+    long n = N.value;
+    double s = 0;
+    for(; e > E.limit; e += E.delta) {
+      s = Modular.addMod(s, montgomery.set(n).mod(e)/(double)n);
+      n += N.delta;
+    }
+    return s;
+  }
+
+  /** {@inheritDoc} */
+  @Override
+  public int compareTo(Summation that) {
+    final int de = this.E.compareTo(that.E);
+    if (de != 0) return de;
+    return this.N.compareTo(that.N);
+  }
+
+  /** {@inheritDoc} */
+  @Override
+  public Summation combine(Summation that) {
+    if (this.N.delta != that.N.delta || this.E.delta != that.E.delta) 
+      throw new IllegalArgumentException(
+          "this.N.delta != that.N.delta || this.E.delta != that.E.delta"
+          + ",\n  this=" + this
+          + ",\n  that=" + that);
+    if (this.E.limit == that.E.value && this.N.limit == that.N.value) {
+      final double v = Modular.addMod(this.value, that.value);
+      final Summation s = new Summation(
+          new ArithmeticProgression(N.symbol, N.value, N.delta, that.N.limit),
+          new ArithmeticProgression(E.symbol, E.value, E.delta, that.E.limit));
+      s.setValue(v);
+      return s;
+    }
+    return null;
+  }
+
+  /** Find the remaining terms. */
+  public <T extends Container<Summation>> List<Summation> remainingTerms(List<T> sorted) {
+    final List<Summation> results = new ArrayList<Summation>();
+    Summation remaining = this;
+
+    if (sorted != null)
+      for(Container<Summation> c : sorted) {
+        final Summation sigma = c.getElement();
+        if (!remaining.contains(sigma))
+          throw new IllegalArgumentException("!remaining.contains(s),"
+              + "\n  remaining = " + remaining
+              + "\n  s         = " + sigma          
+              + "\n  this      = " + this
+              + "\n  sorted    = " + sorted);
+
+        final Summation s = new Summation(sigma.N.limit, N.delta, remaining.N.limit,
+                                          sigma.E.limit, E.delta, remaining.E.limit);
+        if (s.getSteps() > 0)
+          results.add(s);
+        remaining = new Summation(remaining.N.value, N.delta, sigma.N.value,
+                                  remaining.E.value, E.delta, sigma.E.value);
+      }
+
+    if (remaining.getSteps() > 0)
+      results.add(remaining);
+  
+    return results;
+  }
+
+  /** Does this contains that? */
+  public boolean contains(Summation that) {
+    return this.N.contains(that.N) && this.E.contains(that.E);    
+  }
+
+  /** Partition the summation. */
+  public Summation[] partition(final int nParts) {
+    final Summation[] parts = new Summation[nParts];
+    final long steps = (E.limit - E.value)/E.delta + 1;
+
+    long prevN = N.value;
+    long prevE = E.value;
+
+    for(int i = 1; i < parts.length; i++) {
+      final long k = (i * steps)/parts.length;
+
+      final long currN = N.skip(k);
+      final long currE = E.skip(k);
+      parts[i - 1] = new Summation(
+          new ArithmeticProgression(N.symbol, prevN, N.delta, currN),
+          new ArithmeticProgression(E.symbol, prevE, E.delta, currE));
+
+      prevN = currN;
+      prevE = currE;
+    }
+
+    parts[parts.length - 1] = new Summation(
+        new ArithmeticProgression(N.symbol, prevN, N.delta, N.limit),
+        new ArithmeticProgression(E.symbol, prevE, E.delta, E.limit));
+    return parts;
+  }
+}
\ No newline at end of file

Added: hadoop/mapreduce/trunk/src/test/mapred/org/apache/hadoop/examples/pi/math/TestLongLong.java
URL: http://svn.apache.org/viewvc/hadoop/mapreduce/trunk/src/test/mapred/org/apache/hadoop/examples/pi/math/TestLongLong.java?rev=790015&view=auto
==============================================================================
--- hadoop/mapreduce/trunk/src/test/mapred/org/apache/hadoop/examples/pi/math/TestLongLong.java (added)
+++ hadoop/mapreduce/trunk/src/test/mapred/org/apache/hadoop/examples/pi/math/TestLongLong.java Wed Jul  1 01:35:55 2009
@@ -0,0 +1,73 @@
+/**
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements.  See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership.  The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License.  You may obtain a copy of the License at
+ *
+ *     http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package org.apache.hadoop.examples.pi.math;
+
+import java.math.BigInteger;
+import java.util.Random;
+
+public class TestLongLong extends junit.framework.TestCase {
+  static final Random RAN = new Random(); 
+  static final long MASK = (1L << (LongLong.SIZE >> 1)) - 1;
+
+  static long nextPositiveLong() {
+    return RAN.nextLong() & MASK;
+  }
+  
+  static void verifyMultiplication(long a, long b) {
+    final LongLong ll = LongLong.multiplication(new LongLong(), a, b);
+    final BigInteger bi = BigInteger.valueOf(a).multiply(BigInteger.valueOf(b));
+
+    final String s = String.format("\na = %x\nb = %x\nll= " + ll + "\nbi= " + bi.toString(16) + "\n", a, b);
+    //System.out.println(s);
+    assertEquals(s, bi, ll.toBigInteger());
+  }
+
+  public void testMultiplication() {
+    for(int i = 0; i < 100; i++) {
+      final long a = nextPositiveLong();
+      final long b = nextPositiveLong();
+      verifyMultiplication(a, b);
+    }
+    final long max = Long.MAX_VALUE & MASK;
+    verifyMultiplication(max, max);
+  }
+
+  static void verifyRightShift(long a, long b) {
+    final LongLong ll = new LongLong().set(a, b);
+    final BigInteger bi = ll.toBigInteger();
+    
+    for(int i = 0; i < LongLong.SIZE >> 1; i++) {
+      final long result = ll.shiftRight(i) & MASK;
+      final long expected = bi.shiftRight(i).longValue() & MASK;
+      final String s = String.format("\na = %x\nb = %x\nll= " + ll + "\nbi= " + bi.toString(16) + "\n", a, b);
+      assertEquals(s, expected, result);
+    }
+
+    final String s = String.format("\na = %x\nb = %x\nll= " + ll + "\nbi= " + bi.toString(16) + "\n", a, b);
+    //System.out.println(s);
+    assertEquals(s, bi, ll.toBigInteger());
+  }
+
+  public void testRightShift() {
+    for(int i = 0; i < 1000; i++) {
+      final long a = nextPositiveLong();
+      final long b = nextPositiveLong();
+      verifyMultiplication(a, b);
+    }
+  }
+}

Added: hadoop/mapreduce/trunk/src/test/mapred/org/apache/hadoop/examples/pi/math/TestModular.java
URL: http://svn.apache.org/viewvc/hadoop/mapreduce/trunk/src/test/mapred/org/apache/hadoop/examples/pi/math/TestModular.java?rev=790015&view=auto
==============================================================================
--- hadoop/mapreduce/trunk/src/test/mapred/org/apache/hadoop/examples/pi/math/TestModular.java (added)
+++ hadoop/mapreduce/trunk/src/test/mapred/org/apache/hadoop/examples/pi/math/TestModular.java Wed Jul  1 01:35:55 2009
@@ -0,0 +1,337 @@
+/**
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements.  See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership.  The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License.  You may obtain a copy of the License at
+ *
+ *     http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package org.apache.hadoop.examples.pi.math;
+
+import java.math.BigInteger;
+import java.util.Random;
+
+import org.apache.hadoop.examples.pi.Util.Timer;
+
+public class TestModular extends junit.framework.TestCase { 
+  private static final Random RANDOM = new Random(); 
+  private static final BigInteger TWO = BigInteger.valueOf(2);
+
+
+  static final int DIV_VALID_BIT = 32;
+  static final long DIV_LIMIT = 1L << DIV_VALID_BIT; 
+
+  // return r/n for n > r > 0
+  static long div(long sum, long r, long n) {
+    long q = 0;
+    int i = DIV_VALID_BIT - 1;
+    for(r <<= 1; r < n; r <<= 1) i--;
+//System.out.printf("  r=%d, n=%d, q=%d\n", r, n, q);
+    
+    for(; i >= 0 ;) {
+      r -= n;
+      q |= (1L << i);
+      if (r <= 0) break;
+      for(; r < n; r <<= 1) i--;
+//System.out.printf("  r=%d, n=%d, q=%d\n", r, n, q);
+    }
+
+    sum += q;
+    return sum < DIV_LIMIT? sum: sum - DIV_LIMIT;
+  }
+ 
+  public void testDiv() {
+    for(long n = 2; n < 100; n++)
+      for(long r = 1; r < n; r++) {
+        final long a = div(0, r, n);
+        final long b = (long)((r*1.0/n) * (1L << DIV_VALID_BIT));
+        final String s = String.format("r=%d, n=%d, a=%X, b=%X", r, n, a, b);
+        assertEquals(s, b, a);
+      }
+  }
+
+  static long[][][] generateRN(int nsize, int rsize) {
+    final long[][][] rn = new long[nsize][][];
+
+    for(int i = 0; i < rn.length; i++) {
+      rn[i] = new long[rsize + 1][];
+      long n = RANDOM.nextLong() & 0xFFFFFFFFFFFFFFFL; 
+      if (n <= 1) n = 0xFFFFFFFFFFFFFFFL - n;
+      rn[i][0] = new long[]{n};
+      final BigInteger N = BigInteger.valueOf(n); 
+
+      for(int j = 1; j < rn[i].length; j++) {
+        long r = RANDOM.nextLong();
+        if (r < 0) r = -r;
+        if (r >= n) r %= n;
+        final BigInteger R = BigInteger.valueOf(r); 
+        rn[i][j] = new long[]{r, R.multiply(R).mod(N).longValue()};
+      }
+    }
+    return rn;
+  }
+
+  static long square_slow(long z, final long n) {
+    long r = 0;
+    for(long s = z; z > 0; z >>= 1) {
+      if ((((int)z) & 1) == 1) {
+        r += s;
+        if (r >= n) r -= n;
+      }
+
+      s <<= 1;
+      if (s >= n) s -= n;
+    }
+    return r;
+  }
+
+  //0 <= r < n < max/2
+  static long square(long r, final long n, long r2p64) {
+    if (r <= Modular.MAX_SQRT_LONG) {
+      r *= r;
+      if (r >= n) r %= n;
+    } else {
+      final int HALF = (63 - Long.numberOfLeadingZeros(n)) >> 1;
+      final int FULL = HALF << 1;
+      final long ONES = (1 << HALF) - 1; 
+  
+      final long high = r >>> HALF;
+      final long low  = r &= ONES;
+
+      r *= r;
+      if (r >= n) r %= n;
+  
+      if (high != 0) {
+        long s = high * high;
+        if (s >= n) s %= n;
+        for(int i = 0; i < FULL; i++)
+          if ((s <<= 1) >= n) s -= n;
+        
+        if (low == 0)
+          r = s;
+        else {
+          long t = high * low;
+          if (t >= n) t %= n;
+          for(int i = -1; i < HALF; i++)
+            if ((t <<= 1) >= n) t -= n;
+          
+          r += s;
+          if (r >= n) r -= n;
+          r += t;
+          if (r >= n) r -= n;
+        }
+      }
+    }
+    return r;
+  }
+  
+  static void squareBenchmarks() {
+    final Timer t = new Timer(false);
+    t.tick("squareBenchmarks(), MAX_SQRT=" + Modular.MAX_SQRT_LONG);
+
+    final long[][][] rn = generateRN(1000, 1000);
+    t.tick("generateRN");
+
+    for(int i = 0; i < rn.length; i++) {
+      final long n = rn[i][0][0];
+      for(int j = 1; j < rn[i].length; j++) {
+        final long r = rn[i][j][0];
+        final long answer = rn[i][j][1];
+        final long s = square_slow(r, n);
+        if (s != answer)
+          assertEquals("r=" + r + ", n=" + n + ", answer=" + answer + " but s=" + s, answer, s);
+      }
+    }
+    t.tick("square_slow");
+
+    for(int i = 0; i < rn.length; i++) {
+      final long n = rn[i][0][0];
+      long r2p64 = (0x4000000000000000L % n) << 1;
+      if (r2p64 >= n) r2p64 -= n;
+      for(int j = 1; j < rn[i].length; j++) {
+        final long r = rn[i][j][0];
+        final long answer = rn[i][j][1];
+        final long s = square(r, n, r2p64);
+        if (s != answer)
+          assertEquals("r=" + r + ", n=" + n + ", answer=" + answer + " but s=" + s, answer, s);
+      }
+    }
+    t.tick("square");
+
+    for(int i = 0; i < rn.length; i++) {
+      final long n = rn[i][0][0];
+      final BigInteger N = BigInteger.valueOf(n);
+      for(int j = 1; j < rn[i].length; j++) {
+        final long r = rn[i][j][0];
+        final long answer = rn[i][j][1];
+        final BigInteger R = BigInteger.valueOf(r);
+        final long s = R.multiply(R).mod(N).longValue();
+        if (s != answer)
+          assertEquals("r=" + r + ", n=" + n + ", answer=" + answer + " but s=" + s, answer, s);
+      }
+    }
+    t.tick("R.multiply(R).mod(N)");
+
+    for(int i = 0; i < rn.length; i++) {
+      final long n = rn[i][0][0];
+      final BigInteger N = BigInteger.valueOf(n);
+      for(int j = 1; j < rn[i].length; j++) {
+        final long r = rn[i][j][0];
+        final long answer = rn[i][j][1];
+        final BigInteger R = BigInteger.valueOf(r);
+        final long s = R.modPow(TWO, N).longValue();
+        if (s != answer)
+          assertEquals("r=" + r + ", n=" + n + ", answer=" + answer + " but s=" + s, answer, s);
+      }
+    }
+    t.tick("R.modPow(TWO, N)");
+  }
+
+  static long[][][] generateEN(int nsize, int esize) {
+    final long[][][] en = new long[nsize][][];
+
+    for(int i = 0; i < en.length; i++) {
+      en[i] = new long[esize + 1][];
+      long n = (RANDOM.nextLong() & 0xFFFFFFFFFFFFFFFL) | 1L; 
+      if (n == 1) n = 3;
+      en[i][0] = new long[]{n};
+      final BigInteger N = BigInteger.valueOf(n); 
+
+      for(int j = 1; j < en[i].length; j++) {
+        long e = RANDOM.nextLong();
+        if (e < 0) e = -e;
+        final BigInteger E = BigInteger.valueOf(e); 
+        en[i][j] = new long[]{e, TWO.modPow(E, N).longValue()};
+      }
+    }
+    return en;
+  }
+
+  /** Compute $2^e \mod n$ for e > 0, n > 2 */
+  static long modBigInteger(final long e, final long n) {
+    long mask = (e & 0xFFFFFFFF00000000L) == 0 ? 0x00000000FFFFFFFFL
+        : 0xFFFFFFFF00000000L;
+    mask &= (e & 0xFFFF0000FFFF0000L & mask) == 0 ? 0x0000FFFF0000FFFFL
+        : 0xFFFF0000FFFF0000L;
+    mask &= (e & 0xFF00FF00FF00FF00L & mask) == 0 ? 0x00FF00FF00FF00FFL
+        : 0xFF00FF00FF00FF00L;
+    mask &= (e & 0xF0F0F0F0F0F0F0F0L & mask) == 0 ? 0x0F0F0F0F0F0F0F0FL
+        : 0xF0F0F0F0F0F0F0F0L;
+    mask &= (e & 0xCCCCCCCCCCCCCCCCL & mask) == 0 ? 0x3333333333333333L
+        : 0xCCCCCCCCCCCCCCCCL;
+    mask &= (e & 0xAAAAAAAAAAAAAAAAL & mask) == 0 ? 0x5555555555555555L
+        : 0xAAAAAAAAAAAAAAAAL;
+
+    final BigInteger N = BigInteger.valueOf(n);
+    long r = 2;
+    for (mask >>= 1; mask > 0; mask >>= 1) {
+      if (r <= Modular.MAX_SQRT_LONG) {
+        r *= r;
+        if (r >= n) r %= n;
+      } else {
+        final BigInteger R = BigInteger.valueOf(r);
+        r = R.multiply(R).mod(N).longValue();
+      }
+
+      if ((e & mask) != 0) {
+        r <<= 1;
+        if (r >= n) r -= n;
+      }
+    }
+    return r;
+  }
+
+  static class Montgomery2 extends Montgomery {
+    /** Compute 2^y mod N for N odd. */
+    long mod2(final long y) {
+      long r0 = R - N; 
+      long r1 = r0 << 1;
+      if (r1 >= N) r1 -= N;
+      
+      for(long mask = Long.highestOneBit(y); mask > 0; mask >>>= 1) {
+        if ((mask & y) == 0) {
+          r1 = product.m(r0, r1);
+          r0 = product.m(r0, r0);
+        } else {
+          r0 = product.m(r0, r1);
+          r1 = product.m(r1, r1);
+        }
+      }
+      return product.m(r0, 1);
+    }
+  }
+  
+  static void modBenchmarks() {
+    final Timer t = new Timer(false);
+    t.tick("modBenchmarks()");
+
+    final long[][][] en = generateEN(10000, 10);
+    t.tick("generateEN");
+
+    for(int i = 0; i < en.length; i++) {
+      final long n = en[i][0][0];
+      for(int j = 1; j < en[i].length; j++) {
+        final long e = en[i][j][0];
+        final long answer = en[i][j][1];
+        final long s = Modular.mod(e, n);
+        if (s != answer)
+          assertEquals("e=" + e + ", n=" + n + ", answer=" + answer + " but s=" + s, answer, s);
+      }
+    }
+    t.tick("Modular.mod");
+    
+    final Montgomery2 m2 = new Montgomery2();
+    for(int i = 0; i < en.length; i++) {
+      final long n = en[i][0][0];
+      m2.set(n);
+      for(int j = 1; j < en[i].length; j++) {
+        final long e = en[i][j][0];
+        final long answer = en[i][j][1];
+        final long s = m2.mod(e);
+        if (s != answer)
+          assertEquals("e=" + e + ", n=" + n + ", answer=" + answer + " but s=" + s, answer, s);
+      }
+    }
+    t.tick("montgomery.mod");
+
+    for(int i = 0; i < en.length; i++) {
+      final long n = en[i][0][0];
+      m2.set(n);
+      for(int j = 1; j < en[i].length; j++) {
+        final long e = en[i][j][0];
+        final long answer = en[i][j][1];
+        final long s = m2.mod2(e);
+        if (s != answer)
+          assertEquals("e=" + e + ", n=" + n + ", answer=" + answer + " but s=" + s, answer, s);
+      }
+    }
+    t.tick("montgomery.mod2");
+
+    for(int i = 0; i < en.length; i++) {
+      final long n = en[i][0][0];
+      final BigInteger N = BigInteger.valueOf(n); 
+      for(int j = 1; j < en[i].length; j++) {
+        final long e = en[i][j][0];
+        final long answer = en[i][j][1];
+        final long s = TWO.modPow(BigInteger.valueOf(e), N).longValue();
+        if (s != answer)
+          assertEquals("e=" + e + ", n=" + n + ", answer=" + answer + " but s=" + s, answer, s);
+      }
+    }
+    t.tick("BigInteger.modPow(e, n)");
+  }
+
+  public static void main(String[] args) {
+    squareBenchmarks();
+    modBenchmarks();
+  }
+}
\ No newline at end of file

Added: hadoop/mapreduce/trunk/src/test/mapred/org/apache/hadoop/examples/pi/math/TestSummation.java
URL: http://svn.apache.org/viewvc/hadoop/mapreduce/trunk/src/test/mapred/org/apache/hadoop/examples/pi/math/TestSummation.java?rev=790015&view=auto
==============================================================================
--- hadoop/mapreduce/trunk/src/test/mapred/org/apache/hadoop/examples/pi/math/TestSummation.java (added)
+++ hadoop/mapreduce/trunk/src/test/mapred/org/apache/hadoop/examples/pi/math/TestSummation.java Wed Jul  1 01:35:55 2009
@@ -0,0 +1,147 @@
+/**
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements.  See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership.  The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License.  You may obtain a copy of the License at
+ *
+ *     http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package org.apache.hadoop.examples.pi.math;
+
+import java.math.BigInteger;
+import java.util.ArrayList;
+import java.util.Arrays;
+import java.util.Collections;
+import java.util.List;
+import java.util.Random;
+
+import org.apache.hadoop.examples.pi.Container;
+import org.apache.hadoop.examples.pi.Util;
+import org.apache.hadoop.examples.pi.Util.Timer;
+import org.apache.hadoop.examples.pi.math.TestModular.Montgomery2;
+
+public class TestSummation extends junit.framework.TestCase {
+  static final Random RANDOM = new Random();
+  static final BigInteger TWO = BigInteger.valueOf(2);
+
+  private static Summation2 newSummation(final long base, final long range, final long delta) {
+    final ArithmeticProgression N = new ArithmeticProgression('n', base+3, delta, base+3+range);
+    final ArithmeticProgression E = new ArithmeticProgression('e', base+range, -delta, base);
+    return new Summation2(N, E);
+  }
+
+  private static void runTestSubtract(Summation sigma, List<Summation> diff) {
+//    Util.out.println("diff=" + diff);
+    List<Container<Summation>> tmp = new ArrayList<Container<Summation>>(diff.size());
+    for(Summation s : diff)
+      tmp.add(s);
+    final List<Summation> a = sigma.remainingTerms(tmp);
+//    Util.out.println("a   =" + a);
+
+    a.addAll(diff);
+    for(Summation s : a)
+      s.compute();
+
+    final List<Summation> combined = Util.combine(a);
+//    Util.out.println("combined=" + combined);
+    assertEquals(1, combined.size());
+    assertEquals(sigma, combined.get(0));
+  }
+
+  public void testSubtract() {
+    final Summation sigma = newSummation(3, 10000, 20);
+    final int size = 10;
+    final List<Summation> parts = Arrays.asList(sigma.partition(size));
+    Collections.sort(parts);
+
+    runTestSubtract(sigma, new ArrayList<Summation>());
+    runTestSubtract(sigma, parts);
+
+    for(int n = 1; n < size; n++) {
+      for(int j = 0; j < 10; j++) {
+        final List<Summation> diff = new ArrayList<Summation>(parts);
+        for(int i = 0; i < n; i++)
+          diff.remove(RANDOM.nextInt(diff.size()));
+///        Collections.sort(diff);
+        runTestSubtract(sigma, diff);
+      }
+    }
+  }
+  
+  static class Summation2 extends Summation {
+    Summation2(ArithmeticProgression N, ArithmeticProgression E) {
+      super(N, E);
+    }
+
+    
+    final Montgomery2 m2 = new Montgomery2();
+    double compute_montgomery2() {
+      long e = E.value;
+      long n = N.value;
+      double s = 0;
+      for(; e > E.limit; e += E.delta) {
+        m2.set(n);
+        s = Modular.addMod(s, m2.mod2(e)/(double)n);
+        n += N.delta;
+      }
+      return s;
+    }
+
+    double compute_modBigInteger() {
+      long e = E.value;
+      long n = N.value;
+      double s = 0;
+      for(; e > E.limit; e += E.delta) {
+        s = Modular.addMod(s, TestModular.modBigInteger(e, n)/(double)n);
+        n += N.delta;
+      }
+      return s;
+    }
+  
+    double compute_modPow() {
+      long e = E.value;
+      long n = N.value;
+      double s = 0;
+      for(; e > E.limit; e += E.delta) {
+        s = Modular.addMod(s, TWO.modPow(BigInteger.valueOf(e), BigInteger.valueOf(n)).doubleValue()/n);
+        n += N.delta;
+      }
+      return s;
+    }
+  }
+
+  private static void computeBenchmarks(final Summation2 sigma) {
+    final Timer t = new Timer(false);
+    t.tick("sigma=" + sigma);
+    final double value = sigma.compute();
+    t.tick("compute=" + value);
+    assertEquals(value, sigma.compute_modular());
+    t.tick("compute_modular");
+    assertEquals(value, sigma.compute_montgomery());
+    t.tick("compute_montgomery");
+    assertEquals(value, sigma.compute_montgomery2());
+    t.tick("compute_montgomery2");
+
+    assertEquals(value, sigma.compute_modBigInteger());
+    t.tick("compute_modBigInteger");
+    assertEquals(value, sigma.compute_modPow());
+    t.tick("compute_modPow");
+  }
+
+  /** Benchmarks */
+  public static void main(String[] args) {
+    final long delta = 1L << 4;
+    final long range = 1L << 20;
+    for(int i = 20; i < 40; i += 2) 
+      computeBenchmarks(newSummation(1L << i, range, delta));
+  }
+}