commons-dev mailing list archives

Site index · List index
Message view « Date » · « Thread »
Top « Date » · « Thread »
From pste...@apache.org
Subject svn commit: r476930 [14/22] - in /jakarta/commons/proper/math/trunk/src/mantissa: ./ src/ src/org/ src/org/spaceroots/ src/org/spaceroots/mantissa/ src/org/spaceroots/mantissa/algebra/ src/org/spaceroots/mantissa/estimation/ src/org/spaceroots/mantissa...
Date Sun, 19 Nov 2006 21:36:23 GMT
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/GaussLegendreIntegrator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/GaussLegendreIntegrator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/GaussLegendreIntegrator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/GaussLegendreIntegrator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,162 @@
+// 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.spaceroots.mantissa.quadrature.vectorial;
+
+import org.spaceroots.mantissa.functions.vectorial.ComputableFunction;
+import org.spaceroots.mantissa.functions.FunctionException;
+
+/** This class implements a Gauss-Legendre integrator.
+
+ * <p>Gauss-Legendre integrators are efficient integrators that can
+ * accurately integrate functions with few functions evaluations. A
+ * Gauss-Legendre integrator using an n-points quadrature formula can
+ * integrate exactly 2n-1 degree polynoms.</p>
+
+ * <p>These integrators evaluate the function on n carefully chosen
+ * points in each step interval. These points are not evenly
+ * spaced. The function is <emph>never</emph> evaluated at the
+ * boundary points, which means it can be undefined at these
+ * points.</p>
+
+ * @version $Id: GaussLegendreIntegrator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class GaussLegendreIntegrator
+  implements ComputableFunctionIntegrator {
+  /** Build a Gauss-Legendre integrator.
+
+   * <p>A Gauss-Legendre integrator is a formula like:
+   * <pre>
+   *    int (f) from -1 to +1 = Sum (ai * f(xi))
+   * </pre>
+   * </p>
+   *
+   * <p>The coefficients of the formula are computed as follow:
+   * <pre>
+   *   let n be the desired number of points
+   *   the xi are the roots of the degree n Legendre polynomial
+   *   the ai are the integrals int (Li^2) from -1 to +1
+   *   where Li (x) = Prod (x-xk)/(xi-xk) for k != i
+   * </pre>
+   * </p>
+   *
+   * <p>A formula in n points can integrate exactly polynoms of degree
+   * up to 2n-1.</p>
+   *
+   * @param minPoints minimal number of points desired
+   * @param rawStep raw integration step (the precise step will be
+   * adjusted in order to have an integer number of steps in the
+   * integration range).
+   * */
+  public GaussLegendreIntegrator(int minPoints, double rawStep) {
+    if (minPoints <= 2) {
+      weightedRoots = new double[][] {
+        { 1.0, -1.0 / Math.sqrt (3.0) },
+        { 1.0,  1.0 / Math.sqrt (3.0) }
+      };
+    } else if (minPoints <= 3) {
+      weightedRoots = new double[][] {
+        { 5.0 / 9.0, -Math.sqrt (0.6) },
+        { 8.0 / 9.0,             0.0  },
+        { 5.0 / 9.0,  Math.sqrt (0.6) }
+      };
+    } else if (minPoints <= 4) {
+      weightedRoots = new double[][] {
+        { (90.0 - 5.0 * Math.sqrt (30.0)) / 180.0,
+             -Math.sqrt ((15.0 + 2.0 * Math.sqrt (30.0)) / 35.0) },
+        { (90.0 + 5.0 * Math.sqrt (30.0)) / 180.0,
+             -Math.sqrt ((15.0 - 2.0 * Math.sqrt (30.0)) / 35.0) },
+        { (90.0 + 5.0 * Math.sqrt (30.0)) / 180.0,
+              Math.sqrt ((15.0 - 2.0 * Math.sqrt (30.0)) / 35.0) },
+        { (90.0 - 5.0 * Math.sqrt (30.0)) / 180.0,
+              Math.sqrt ((15.0 + 2.0 * Math.sqrt (30.0)) / 35.0) }
+      };
+    } else {
+      weightedRoots = new double[][] {
+        { (322.0 - 13.0 * Math.sqrt (70.0)) / 900.0,
+             -Math.sqrt ((35.0 + 2.0 * Math.sqrt (70.0)) / 63.0) },
+        { (322.0 + 13.0 * Math.sqrt (70.0)) / 900.0,
+             -Math.sqrt ((35.0 - 2.0 * Math.sqrt (70.0)) / 63.0) },
+        { 128.0 / 225.0,
+              0.0 },
+        { (322.0 + 13.0 * Math.sqrt (70.0)) / 900.0,
+              Math.sqrt ((35.0 - 2.0 * Math.sqrt (70.0)) / 63.0) },
+        { (322.0 - 13.0 * Math.sqrt (70.0)) / 900.0,
+              Math.sqrt ((35.0 + 2.0 * Math.sqrt (70.0)) / 63.0) }
+      };
+    }
+
+    this.rawStep = rawStep;
+
+  }
+
+  /** Get the number of functions evaluation per step.
+   * @return number of function evaluation per step
+   */
+  public int getEvaluationsPerStep() {
+    return weightedRoots.length;
+  }
+
+  public double[] integrate(ComputableFunction f, double a, double b)
+    throws FunctionException {
+
+    // swap the integration bounds if they are not in ascending order
+    if (b < a) {
+      double tmp = b;
+      b          = a;
+      a          = tmp;
+    }
+
+    // adjust the integration step according to the bounds
+    long   n     = Math.round(0.5 + (b - a) / rawStep);
+    double step  = (b - a) / n;
+
+    // integrate over all elementary steps
+    double halfStep = step / 2.0;
+    double midPoint = a + halfStep;
+
+    double[] sum = new double[f.getDimension()];
+    for (int k = 0; k < sum.length; ++k) {
+      sum[k] = 0.0;
+    }
+
+    for (long i = 0; i < n; ++i) {
+      for (int j = 0; j < weightedRoots.length; ++j) {
+        double[] value = f.valueAt(midPoint + halfStep * weightedRoots[j][1]);
+        for (int k = 0; k < sum.length; ++k) {
+          sum[k] += weightedRoots[j][0] * value[k];
+        }
+      }
+      midPoint += step;
+    }
+
+    for (int k = 0; k < sum.length; ++k) {
+      sum [k] *= halfStep;
+    }
+
+    return sum;
+
+  }
+
+  double[][] weightedRoots;
+
+  double rawStep;
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/GaussLegendreIntegrator.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/RiemannIntegrator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/RiemannIntegrator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/RiemannIntegrator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/RiemannIntegrator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,63 @@
+// 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.spaceroots.mantissa.quadrature.vectorial;
+
+import org.spaceroots.mantissa.functions.FunctionException;
+import org.spaceroots.mantissa.functions.ExhaustedSampleException;
+import org.spaceroots.mantissa.functions.vectorial.SampledFunctionIterator;
+
+/** This class implements a Riemann integrator.
+
+ * <p>A Riemann integrator is a very simple one that assumes the
+ * function is constant over the integration step. Since it is very
+ * simple, this algorithm needs very small steps to achieve high
+ * accuracy, and small steps lead to numerical errors and
+ * instabilities.</p>
+
+ * <p>This algorithm is almost never used and has been included in
+ * this package only as a simple template for more useful
+ * integrators.</p>
+
+ * @see TrapezoidIntegrator
+
+ * @version $Id: RiemannIntegrator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class RiemannIntegrator
+  implements SampledFunctionIntegrator {
+
+  public double[] integrate(SampledFunctionIterator iter)
+    throws ExhaustedSampleException, FunctionException {
+
+    RiemannIntegratorSampler sampler = new RiemannIntegratorSampler(iter);
+    double[] sum = null;
+
+    try {
+      while (true) {
+        sum = sampler.nextSamplePoint().getY();
+      }
+    } catch(ExhaustedSampleException e) {
+    }
+
+    return sum;
+
+  }
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/RiemannIntegrator.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/RiemannIntegratorSampler.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/RiemannIntegratorSampler.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/RiemannIntegratorSampler.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/RiemannIntegratorSampler.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,101 @@
+// 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.spaceroots.mantissa.quadrature.vectorial;
+
+import org.spaceroots.mantissa.functions.vectorial.*;
+import org.spaceroots.mantissa.functions.FunctionException;
+import org.spaceroots.mantissa.functions.ExhaustedSampleException;
+
+/** This class implements a Riemann integrator as a sample.
+
+ * <p>A Riemann integrator is a very simple one that assumes the
+ * function is constant over the integration step. Since it is very
+ * simple, this algorithm needs very small steps to achieve high
+ * accuracy, and small steps lead to numerical errors and
+ * instabilities.</p>
+
+ * <p>This algorithm is almost never used and has been included in
+ * this package only as a simple template for more useful
+ * integrators.</p>
+
+ * @see RiemannIntegrator
+
+ * @version $Id: RiemannIntegratorSampler.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class RiemannIntegratorSampler
+  implements SampledFunctionIterator {
+
+  /** Underlying sample iterator. */
+  private SampledFunctionIterator iter;
+
+  /** Current point. */
+  private VectorialValuedPair current;
+
+  /** Current running sum. */
+  private double[] sum;
+
+  /** Constructor.
+   * Build an integrator from an underlying sample iterator.
+   * @param iter iterator over the base function
+   */
+  public RiemannIntegratorSampler(SampledFunctionIterator iter)
+    throws ExhaustedSampleException, FunctionException {
+
+    this.iter = iter;
+
+    // get the first point
+    current = iter.nextSamplePoint();
+
+    // initialize the sum
+    sum = new double[iter.getDimension()];
+    for (int i = 0; i < sum.length; ++i) {
+      sum[i] = 0.0;
+    }
+
+  }
+
+  public boolean hasNext() {
+    return iter.hasNext();
+  }
+
+  public int getDimension() {
+    return iter.getDimension();
+  }
+
+  public VectorialValuedPair nextSamplePoint()
+    throws ExhaustedSampleException, FunctionException {
+
+    // performs one step of a Riemann scheme
+    VectorialValuedPair previous = current;
+    current = iter.nextSamplePoint();
+    double step =  (current.getX() - previous.getX());
+    double[] pY = previous.getY();
+    for (int i = 0; i < sum.length; ++i) {
+      sum[i] += step * pY[i];
+    }
+
+    double[] values = new double[sum.length];
+    System.arraycopy(sum, 0, values, 0, sum.length);
+    return new VectorialValuedPair (current.getX(), values);
+
+  }
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/RiemannIntegratorSampler.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/SampledFunctionIntegrator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/SampledFunctionIntegrator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/SampledFunctionIntegrator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/SampledFunctionIntegrator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,49 @@
+// 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.spaceroots.mantissa.quadrature.vectorial;
+
+import org.spaceroots.mantissa.functions.vectorial.SampledFunctionIterator;
+import org.spaceroots.mantissa.functions.ExhaustedSampleException;
+import org.spaceroots.mantissa.functions.FunctionException;
+
+/** This interface represents an integrator for vectorial samples.
+
+ * <p>The classes which are devoted to integrate vectorial samples
+ * should implement this interface.</p>
+
+ * @see org.spaceroots.mantissa.functions.vectorial.SampledFunctionIterator
+ * @see ComputableFunctionIntegrator
+
+ * @version $Id: SampledFunctionIntegrator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public interface SampledFunctionIntegrator {
+
+  /** Integrate a sample over its overall range
+   * @param iter iterator over the sample to integrate
+   * @return value of the integral over the sample range
+   * @exception ExhaustedSampleException if the sample does not have
+   * enough points for the integration scheme
+   * @exception FunctionException if the underlying sampled function throws one
+   */
+  public double[] integrate(SampledFunctionIterator iter)
+    throws ExhaustedSampleException, FunctionException;
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/SampledFunctionIntegrator.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/TrapezoidIntegrator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/TrapezoidIntegrator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/TrapezoidIntegrator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/TrapezoidIntegrator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,54 @@
+// 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.spaceroots.mantissa.quadrature.vectorial;
+
+import org.spaceroots.mantissa.functions.FunctionException;
+import org.spaceroots.mantissa.functions.ExhaustedSampleException;
+import org.spaceroots.mantissa.functions.vectorial.SampledFunctionIterator;
+
+/** This class implements a trapezoid integrator.
+
+ * <p>A trapezoid integrator is a very simple one that assumes the
+ * function is linear over the integration step.</p>
+
+ * @version $Id: TrapezoidIntegrator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class TrapezoidIntegrator
+  implements SampledFunctionIntegrator {
+  public double[] integrate(SampledFunctionIterator iter)
+    throws ExhaustedSampleException, FunctionException {
+
+    TrapezoidIntegratorSampler sampler =
+      new TrapezoidIntegratorSampler(iter);
+    double[] sum = null;
+
+    try {
+      while (true) {
+        sum = sampler.nextSamplePoint().getY();
+      }
+    } catch(ExhaustedSampleException e) {
+    }
+
+    return sum;
+
+  }
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/TrapezoidIntegrator.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/TrapezoidIntegratorSampler.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/TrapezoidIntegratorSampler.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/TrapezoidIntegratorSampler.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/TrapezoidIntegratorSampler.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,103 @@
+// 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.spaceroots.mantissa.quadrature.vectorial;
+
+import org.spaceroots.mantissa.functions.vectorial.*;
+import org.spaceroots.mantissa.functions.FunctionException;
+import org.spaceroots.mantissa.functions.ExhaustedSampleException;
+
+/** This class implements a trapezoid integrator as a sample.
+
+ * <p>A trapezoid integrator is a very simple one that assumes the
+ * function is constant over the integration step. Since it is very
+ * simple, this algorithm needs very small steps to achieve high
+ * accuracy, and small steps lead to numerical errors and
+ * instabilities.</p>
+
+ * <p>This algorithm is almost never used and has been included in
+ * this package only as a simple template for more useful
+ * integrators.</p>
+
+ * @see TrapezoidIntegrator
+
+ * @version $Id: TrapezoidIntegratorSampler.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class TrapezoidIntegratorSampler
+  implements SampledFunctionIterator {
+
+  /** Underlying sample iterator. */
+  private SampledFunctionIterator iter;
+
+  /** Current point. */
+  private VectorialValuedPair current;
+
+  /** Current running sum. */
+  private double[] sum;
+
+  /** Constructor.
+   * Build an integrator from an underlying sample iterator.
+   * @param iter iterator over the base function
+   */
+  public TrapezoidIntegratorSampler(SampledFunctionIterator iter)
+    throws ExhaustedSampleException, FunctionException {
+
+    this.iter = iter;
+
+    // get the first point
+    current = iter.nextSamplePoint();
+
+    // initialize the sum
+    sum = new double[iter.getDimension()];
+    for (int i = 0; i < sum.length; ++i) {
+      sum[i] = 0.0;
+    }
+
+  }
+
+  public boolean hasNext() {
+    return iter.hasNext();
+  }
+
+  public int getDimension() {
+    return iter.getDimension();
+  }
+
+  public VectorialValuedPair nextSamplePoint()
+    throws ExhaustedSampleException, FunctionException {
+
+    // performs one step of a trapezoid scheme
+    VectorialValuedPair previous = current;
+    current = iter.nextSamplePoint();
+
+    double halfDx = 0.5 * (current.getX() - previous.getX());
+    double[] pY = previous.getY();
+    double[] cY = current.getY();
+    for (int i = 0; i < sum.length; ++i) {
+      sum[i] += halfDx * (pY[i] + cY[i]);
+    }
+
+    double[] values = new double[sum.length];
+    System.arraycopy(sum, 0, values, 0, sum.length);
+    return new VectorialValuedPair (current.getX(), values);
+
+  }
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/TrapezoidIntegratorSampler.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/CorrelatedRandomVectorGenerator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/CorrelatedRandomVectorGenerator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/CorrelatedRandomVectorGenerator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/CorrelatedRandomVectorGenerator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,287 @@
+// 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.spaceroots.mantissa.random;
+
+import org.spaceroots.mantissa.MantissaException;
+import org.spaceroots.mantissa.linalg.Matrix;
+import org.spaceroots.mantissa.linalg.GeneralMatrix;
+import org.spaceroots.mantissa.linalg.SymetricalMatrix;
+
+import java.io.Serializable;
+
+/** This class allows to generate random vectors with correlated components.
+
+ * <p>Random vectors with correlated components are built by combining
+ * the uncorrelated components of another random vector in such a way
+ * the resulting correlations are the ones specified by a positive
+ * definite covariance matrix.</p>
+
+ * <p>Sometimes, the covariance matrix for a given simulation is not
+ * strictly positive definite. This means that the correlations are
+ * not all independant from each other. In this case, however, the non
+ * strictly positive elements found during the Cholesky decomposition
+ * of the covariance matrix should not be negative either, they
+ * should be null. This implies that rather than computing <code>C =
+ * L.Lt</code> where <code>C</code> is the covariance matrix and
+ * <code>L</code> is a lower-triangular matrix, we compute <code>C =
+ * B.Bt</code> where <code>B</code> is a rectangular matrix having
+ * more rows than columns. The number of columns of <code>B</code> is
+ * the rank of the covariance matrix, and it is the dimension of the
+ * uncorrelated random vector that is needed to compute the component
+ * of the correlated vector. This class does handle this situation
+ * automatically.</p>
+
+ * @version $Id: CorrelatedRandomVectorGenerator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class CorrelatedRandomVectorGenerator
+  implements Serializable, RandomVectorGenerator {
+
+  /** Simple constructor.
+   * <p>Build a correlated random vector generator from its mean
+   * vector and covariance matrix.</p>
+   * @param mean expected mean values for all components
+   * @param covariance covariance matrix
+   * @param generator underlying generator for uncorrelated normalized
+   * components
+   * @exception IllegalArgumentException if there is a dimension
+   * mismatch between the mean vector and the covariance matrix
+   * @exception NotPositiveDefiniteMatrixException if the
+   * covariance matrix is not strictly positive definite
+   */
+  public CorrelatedRandomVectorGenerator(double[] mean,
+                                         SymetricalMatrix covariance,
+                                         NormalizedRandomGenerator generator)
+    throws NotPositiveDefiniteMatrixException {
+
+    int order = covariance.getRows();
+    if (mean.length != order) {
+      String message =
+        MantissaException.translate("dimension mismatch {0} != {1}",
+                                    new String[] {
+                                      Integer.toString(mean.length),
+                                      Integer.toString(order)
+                                    });
+      throw new IllegalArgumentException(message);
+    }
+    this.mean = mean;
+
+    factorize(covariance);
+
+    this.generator = generator;
+    normalized = new double[rank];
+    correlated = new double[order];
+
+  }
+
+  /** Simple constructor.
+   * <p>Build a null mean random correlated vector generator from its
+   * covariance matrix.</p>
+   * @param covariance covariance matrix
+   * @param generator underlying generator for uncorrelated normalized
+   * components
+   * @exception NotPositiveDefiniteMatrixException if the
+   * covariance matrix is not strictly positive definite
+   */
+  public CorrelatedRandomVectorGenerator(SymetricalMatrix covariance,
+                                         NormalizedRandomGenerator generator)
+    throws NotPositiveDefiniteMatrixException {
+
+    int order = covariance.getRows();
+    mean = new double[order];
+    for (int i = 0; i < order; ++i) {
+      mean[i] = 0;
+    }
+
+    factorize(covariance);
+
+    this.generator = generator;
+    normalized = new double[rank];
+    correlated = new double[order];
+
+  }
+
+  /** Get the root of the covariance matrix.
+   * The root is the matrix <code>B</code> such that <code>B.Bt</code>
+   * is equal to the covariance matrix
+   * @return root of the square matrix
+   */
+  public Matrix getRootMatrix() {
+    return root;
+  }
+
+  /** Get the underlying normalized components generator.
+   * @return underlying uncorrelated components generator
+   */
+  public NormalizedRandomGenerator getGenerator() {
+    return generator;
+  }
+
+  /** Get the rank of the covariance matrix.
+   * The rank is the number of independant rows in the covariance
+   * matrix, it is also the number of columns of the rectangular
+   * matrix of the factorization.
+   * @return rank of the square matrix.
+   */
+  public int getRank() {
+    return rank;
+  }
+
+  /** Factorize the original square matrix.
+   * @param covariance covariance matrix
+   * @exception NotPositiveDefiniteMatrixException if the
+   * covariance matrix is not strictly positive definite
+   */
+  private void factorize(SymetricalMatrix covariance)
+  throws NotPositiveDefiniteMatrixException {
+
+    int order = covariance.getRows();
+    SymetricalMatrix c = (SymetricalMatrix) covariance.duplicate();
+    GeneralMatrix    b = new GeneralMatrix(order, order);
+
+    int[] swap  = new int[order];
+    int[] index = new int[order];
+    for (int i = 0; i < order; ++i) {
+      index[i] = i;
+    }
+
+    rank = 0;
+    for (boolean loop = true; loop;) {
+
+      // find maximal diagonal element
+      swap[rank] = rank;
+      for (int i = rank + 1; i < order; ++i) {
+        if (c.getElement(index[i], index[i])
+            > c.getElement(index[swap[i]], index[swap[i]])) {
+          swap[rank] = i;
+        }
+      }
+
+
+      // swap elements
+      if (swap[rank] != rank) {
+        int tmp = index[rank];
+        index[rank] = index[swap[rank]];
+        index[swap[rank]] = tmp;
+      }
+
+      // check diagonal element
+      if (c.getElement(index[rank], index[rank]) < 1.0e-12) {
+
+        if (rank == 0) {
+          throw new NotPositiveDefiniteMatrixException();
+        }
+
+        // check remaining diagonal elements
+        for (int i = rank; i < order; ++i) {
+          if (c.getElement(index[rank], index[rank]) < -1.0e-12) {
+            // there is at least one sufficiently negative diagonal element,
+            // the covariance matrix is wrong
+            throw new NotPositiveDefiniteMatrixException();
+          }
+        }
+
+        // all remaining diagonal elements are close to zero,
+        // we consider we have found the rank of the covariance matrix
+        ++rank;
+        loop = false;
+
+      } else {
+
+        // transform the matrix
+        double sqrt = Math.sqrt(c.getElement(index[rank], index[rank]));
+        b.setElement(rank, rank, sqrt);
+        double inverse = 1 / sqrt;
+        for (int i = rank + 1; i < order; ++i) {
+          double e = inverse * c.getElement(index[i], index[rank]);
+          b.setElement(i, rank, e);
+          c.setElement(index[i], index[i],
+                          c.getElement(index[i], index[i]) - e * e);
+          for (int j = rank + 1; j < i; ++j) {
+            double f = b.getElement(j, rank);
+            c.setElementAndSymetricalElement(index[i], index[j],
+                                             c.getElement(index[i], index[j])
+                                             - e * f);
+          }
+        }
+
+        // prepare next iteration
+        loop = ++rank < order;
+
+      }
+
+    }
+
+    // build the root matrix
+    root = new GeneralMatrix(order, rank);
+    for (int i = 0; i < order; ++i) {
+      for (int j = 0; j < rank; ++j) {
+        root.setElement(swap[i], j, b.getElement(i, j));
+      }
+    }
+
+  }
+
+  /** Generate a correlated random vector.
+   * @return a random vector as an array of double. The generator
+   * <em>will</em> reuse the same array for each call, in order to
+   * save the allocation time, so the user should keep a copy by
+   * himself if he needs so.
+   */
+  public double[] nextVector() {
+
+    // generate uncorrelated vector
+    for (int i = 0; i < rank; ++i) {
+      normalized[i] = generator.nextDouble();
+    }
+
+    // compute correlated vector
+    for (int i = 0; i < correlated.length; ++i) {
+      correlated[i] = mean[i];
+      for (int j = 0; j < rank; ++j) {
+        correlated[i] += root.getElement(i, j) * normalized[j];
+      }
+    }
+
+    return correlated;
+
+  }
+
+  /** Mean vector. */
+  private double[] mean;
+
+  /** Permutated Cholesky root of the covariance matrix. */
+  private Matrix root;
+
+  /** Rank of the covariance matrix. */
+  private int rank;
+
+  /** Underlying generator. */
+  NormalizedRandomGenerator generator;
+
+  /** Storage for the normalized vector. */
+  private double[] normalized;
+
+  /** Storage for the random vector. */
+  private double[] correlated;
+
+  private static final long serialVersionUID = -4754497552287369719L;
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/CorrelatedRandomVectorGenerator.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/GaussianRandomGenerator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/GaussianRandomGenerator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/GaussianRandomGenerator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/GaussianRandomGenerator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,67 @@
+// 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.spaceroots.mantissa.random;
+
+import java.util.Random;
+
+/** This class is a gaussian normalized random generator
+ * for scalars.
+
+ * <p>This class is a simple interface adaptor around the {@link
+ * java.util.Random#nextGaussian nextGaussian} method.</p>
+
+ * @version $Id: GaussianRandomGenerator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class GaussianRandomGenerator
+  implements NormalizedRandomGenerator {
+
+  /** Underlying generator. */
+  Random generator;
+
+  /** Create a new generator.
+   * The seed of the generator is related to the current time.
+   */
+  public GaussianRandomGenerator() {
+    generator = new Random();
+  }
+
+  /** Creates a new random number generator using a single int seed.
+   * @param seed the initial seed (32 bits integer)
+   */
+  public GaussianRandomGenerator(int seed) {
+    generator = new Random(seed);
+  }
+
+  /** Create a new generator initialized with a single long seed.
+   * @param seed seed for the generator (64 bits integer)
+   */
+  public GaussianRandomGenerator(long seed) {
+    generator = new Random(seed);
+  }
+
+  /** Generate a random scalar with null mean and unit standard deviation.
+   * @return a random scalar with null mean and unit standard deviation
+   */
+  public double nextDouble() {
+    return generator.nextGaussian();
+  }
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/GaussianRandomGenerator.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/NormalizedRandomGenerator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/NormalizedRandomGenerator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/NormalizedRandomGenerator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/NormalizedRandomGenerator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,38 @@
+// 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.spaceroots.mantissa.random;
+
+/** This interface represent a normalized random generator for
+ * scalars.
+ * Normalized generator should provide null mean and unit standard
+ * deviation scalars.
+ * @version $Id: NormalizedRandomGenerator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+ */
+public interface NormalizedRandomGenerator {
+
+  /** Generate a random scalar with null mean and unit standard deviation.
+   * <p>This method does <strong>not</strong> specify the shape of the
+   * distribution, it is the implementing class that provides it. The
+   * only contract here is to generate numbers with null mean and unit
+   * standard deviation.</p>
+   * @return a random scalar
+   */
+  public double nextDouble();
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/NormalizedRandomGenerator.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/NotPositiveDefiniteMatrixException.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/NotPositiveDefiniteMatrixException.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/NotPositiveDefiniteMatrixException.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/NotPositiveDefiniteMatrixException.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,50 @@
+// 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.spaceroots.mantissa.random;
+
+import org.spaceroots.mantissa.MantissaException;
+
+/** This class represents exceptions thrown by the correlated random
+ * vector generator.
+
+ * @version $Id: NotPositiveDefiniteMatrixException.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class NotPositiveDefiniteMatrixException
+  extends MantissaException {
+
+  /** Simple constructor.
+   * build an exception with a default message.
+   */
+  public NotPositiveDefiniteMatrixException() {
+    super("not positive definite matrix");
+  }
+
+  /** Simple constructor.
+   * build an exception with the specified message.
+   * @param message message to use to build the exception
+   */
+  public NotPositiveDefiniteMatrixException(String message) {
+    super(message);
+  }
+
+  private static final long serialVersionUID = -6801349873804445905L;
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/NotPositiveDefiniteMatrixException.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/RandomVectorGenerator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/RandomVectorGenerator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/RandomVectorGenerator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/RandomVectorGenerator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,37 @@
+// 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.spaceroots.mantissa.random;
+
+/** This interface represent a random generator for whole vectors.
+
+ * @version $Id: RandomVectorGenerator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public interface RandomVectorGenerator {
+
+  /** Generate a random vector.
+   * @return a random vector as an array of double. The generator
+   * <em>will</em> reuse the same array for each call, in order to
+   * save the allocation time, so the user should keep a copy by
+   * himself if he needs so.
+   */
+  public double[] nextVector();
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/RandomVectorGenerator.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/ScalarSampleStatistics.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/ScalarSampleStatistics.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/ScalarSampleStatistics.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/ScalarSampleStatistics.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,161 @@
+// 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.spaceroots.mantissa.random;
+
+/** This class compute basic statistics on a scalar sample.
+ * @version $Id: ScalarSampleStatistics.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+ */
+public class ScalarSampleStatistics {
+
+  /** Number of sample points. */
+  private int n;
+
+  /** Minimal value in the sample. */
+  private double min;
+
+  /** Maximal value in the sample. */
+  private double max;
+
+  /** Sum of the sample values. */
+  private double sum;
+
+  /** Sum of the squares of the sample values. */
+  private double sum2;
+
+  /** Simple constructor.
+   * Build a new empty instance
+   */
+  public ScalarSampleStatistics() {
+    n    = 0;
+    min  = Double.NaN;
+    max  = min;
+    sum  = 0;
+    sum2 = 0;
+  }
+
+  /** Add one point to the instance.
+   * @param x value of the sample point
+   */
+  public void add(double x) {
+
+    if (n++ == 0) {
+      min  = x;
+      max  = x;
+      sum  = x;
+      sum2 = x * x;
+    } else {
+
+      if (x < min) {
+        min = x;
+      } else if (x > max) {
+        max = x;
+      }
+
+      sum  += x;
+      sum2 += x * x;
+
+    }
+
+  }
+
+  /** Add all points of an array to the instance.
+   * @param points array of points
+   */
+  public void add(double[] points) {
+    for (int i = 0; i < points.length; ++i) {
+      add(points[i]);
+    }
+  }
+
+  /** Add all the points of another sample to the instance.
+   * @param s sample to add
+   */
+  public void add(ScalarSampleStatistics s) {
+
+    if (s.n == 0) {
+      // nothing to add
+      return;
+    }
+
+    if (n == 0) {
+      n    = s.n;
+      min  = s.min;
+      max  = s.max;
+      sum  = s.sum;
+      sum2 = s.sum2;
+    } else {
+
+      n += s.n;
+
+      if (s.min < min) {
+        min = s.min;
+      } else if (s.max > max) {
+        max = s.max;
+      }
+
+      sum  += s.sum;
+      sum2 += s.sum2;
+
+    }
+
+  }
+
+  /** Get the number of points in the sample.
+   * @return number of points in the sample
+   */
+  public int size() {
+    return n;
+  }
+
+  /** Get the minimal value in the sample.
+   * @return minimal value in the sample
+   */
+  public double getMin() {
+    return min;
+  }
+
+  /** Get the maximal value in the sample.
+   * @return maximal value in the sample
+   */
+  public double getMax() {
+    return max;
+  }
+
+  /** Get the mean value of the sample.
+   * @return mean value of the sample
+   */
+  public double getMean() {
+    return (n == 0) ? 0 : (sum / n);
+  }
+
+  /** Get the standard deviation of the underlying probability law.
+   * This method estimate the standard deviation considering that the
+   * data available are only a <em>sample</em> of all possible
+   * values. This value is often called the sample standard deviation
+   * (as opposed to the population standard deviation).
+   * @return standard deviation of the underlying probability law
+   */
+  public double getStandardDeviation() {
+    if (n < 2) {
+      return 0;
+    }
+    return Math.sqrt((n * sum2 - sum * sum) / (n * (n - 1)));
+  }
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/ScalarSampleStatistics.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/UncorrelatedRandomVectorGenerator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/UncorrelatedRandomVectorGenerator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/UncorrelatedRandomVectorGenerator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/UncorrelatedRandomVectorGenerator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,116 @@
+// 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.spaceroots.mantissa.random;
+
+import java.io.Serializable;
+
+/** This class allows to generate random vectors with uncorrelated components.
+
+ * @version $Id: UncorrelatedRandomVectorGenerator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class UncorrelatedRandomVectorGenerator
+  implements Serializable, RandomVectorGenerator {
+
+  /** Simple constructor.
+   * <p>Build an uncorrelated random vector generator from its mean
+   * and standard deviation vectors.</p>
+   * @param mean expected mean values for all components
+   * @param standardDeviation standard deviation for all components
+   * @param generator underlying generator for uncorrelated normalized
+   * components
+   * @exception IllegalArgumentException if there is a dimension
+   * mismatch between the mean and standard deviation vectors
+   */
+  public UncorrelatedRandomVectorGenerator(double[] mean,
+                                           double[] standardDeviation,
+                                           NormalizedRandomGenerator generator) {
+
+    if (mean.length != standardDeviation.length) {
+      throw new IllegalArgumentException("dimension mismatch");
+    }
+    this.mean              = mean;
+    this.standardDeviation = standardDeviation;
+
+    this.generator = generator;
+    random = new double[mean.length];
+
+  }
+
+  /** Simple constructor.
+   * <p>Build a null mean random and unit standard deviation
+   * uncorrelated vector generator</p>
+   * @param dimension dimension of the vectors to generate
+   * @param generator underlying generator for uncorrelated normalized
+   * components
+   */
+  public UncorrelatedRandomVectorGenerator(int dimension,
+                                           NormalizedRandomGenerator generator) {
+
+    mean              = new double[dimension];
+    standardDeviation = new double[dimension];
+    for (int i = 0; i < dimension; ++i) {
+      mean[i]              = 0;
+      standardDeviation[i] = 1;
+    }
+
+    this.generator = generator;
+    random = new double[dimension];
+
+  }
+
+  /** Get the underlying normalized components generator.
+   * @return underlying uncorrelated components generator
+   */
+  public NormalizedRandomGenerator getGenerator() {
+    return generator;
+  }
+
+  /** Generate a correlated random vector.
+   * @return a random vector as an array of double. The generator
+   * <em>will</em> reuse the same array for each call, in order to
+   * save the allocation time, so the user should keep a copy by
+   * himself if he needs so.
+   */
+  public double[] nextVector() {
+
+    for (int i = 0; i < random.length; ++i) {
+      random[i] = mean[i] + standardDeviation[i] * generator.nextDouble();
+    }
+
+    return random;
+
+  }
+
+  /** Mean vector. */
+  private double[] mean;
+
+  /** Standard deviation vector. */
+  private double[] standardDeviation;
+
+  /** Underlying scalar generator. */
+  NormalizedRandomGenerator generator;
+
+  /** Storage for the random vector. */
+  private double[] random;
+
+  private static final long serialVersionUID = -3323293740860311151L;
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/UncorrelatedRandomVectorGenerator.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/UniformRandomGenerator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/UniformRandomGenerator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/UniformRandomGenerator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/UniformRandomGenerator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,74 @@
+// 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.spaceroots.mantissa.random;
+
+import java.util.Random;
+
+/** This class implements a normalized uniform random generator.
+
+ * <p>Since this is a normalized random generator, it has a null mean
+ * and a unit standard deviation. Being also a uniform
+ * generator, it produces numbers in the range [-sqrt(3) ;
+ * sqrt(3)].</p>
+
+ * @version $Id: UniformRandomGenerator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class UniformRandomGenerator
+  implements NormalizedRandomGenerator {
+
+  private static final double SQRT3 = Math.sqrt(3.0);
+
+  private static final double TWOSQRT3 = 2.0 * Math.sqrt(3.0);
+
+  /** Underlying generator. */
+  Random generator;
+
+  /** Create a new generator.
+   * The seed of the generator is related to the current time.
+   */
+  public UniformRandomGenerator() {
+    generator = new Random();
+  }
+
+  /** Creates a new random number generator using a single int seed.
+   * @param seed the initial seed (32 bits integer)
+   */
+  public UniformRandomGenerator(int seed) {
+    generator = new Random(seed);
+  }
+
+  /** Create a new generator initialized with a single long seed.
+   * @param seed seed for the generator (64 bits integer)
+   */
+  public UniformRandomGenerator(long seed) {
+    generator = new Random(seed);
+  }
+
+  /** Generate a random scalar with null mean and unit standard deviation.
+   * <p>The number generated is uniformly distributed between -sqrt(3)
+   * and sqrt(3).</p>
+   * @return a random scalar with null mean and unit standard deviation
+   */
+  public double nextDouble() {
+    return TWOSQRT3 * generator.nextDouble() - SQRT3;
+  }
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/UniformRandomGenerator.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/VectorialSampleStatistics.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/VectorialSampleStatistics.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/VectorialSampleStatistics.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/VectorialSampleStatistics.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,313 @@
+// 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.spaceroots.mantissa.random;
+
+import org.spaceroots.mantissa.linalg.SymetricalMatrix;
+
+import java.util.Arrays;
+
+/** This class compute basic statistics on a scalar sample.
+ * @version $Id: VectorialSampleStatistics.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+ */
+public class VectorialSampleStatistics {
+
+  /** Dimension of the vectors to handle. */
+  private int dimension;
+
+  /** Number of sample points. */
+  private int n;
+
+  /** Indices of the minimal values occurrence in the sample. */
+  private int[] minIndices;
+
+  /** Minimal value in the sample. */
+  private double[] min;
+
+  /** Maximal value in the sample. */
+  private double[] max;
+
+  /** Indices of the maximal values occurrence in the sample. */
+  private int[] maxIndices;
+
+  /** Sum of the sample values. */
+  private double[] sum;
+
+  /** Sum of the squares of the sample values. */
+  private double[] sum2;
+
+  /** Simple constructor.
+   * Build a new empty instance
+   */
+  public VectorialSampleStatistics() {
+    dimension  = -1;
+    n          = 0;
+    min        = null;
+    minIndices = null;
+    max        = null;
+    maxIndices = null;
+    sum        = null;
+    sum2       = null;
+  }
+
+  /** Allocate all the arrays. */
+  private void allocate() {
+    min        = new double[dimension];
+    minIndices = new int[dimension];
+    max        = new double[dimension];
+    maxIndices = new int[dimension];
+    sum        = new double[dimension];
+    sum2       = new double[dimension * (dimension + 1) / 2];
+  }
+
+  /** Add one point to the instance.
+   * @param x value of the sample point
+   * @exception IllegalArgumentException if there is a dimension
+   * mismatch between this point and the ones already added (this
+   * cannot happen when the instance is empty)
+   */
+  public void add(double[] x) {
+
+    if (n == 0) {
+
+      dimension = x.length;
+      allocate();
+
+      Arrays.fill(minIndices, 0);
+      Arrays.fill(maxIndices, 0);
+      System.arraycopy(x, 0, min, 0, dimension);
+      System.arraycopy(x, 0, max, 0, dimension);
+      System.arraycopy(x, 0, sum, 0, dimension);
+
+      int k = 0;
+      for (int i = 0; i < dimension; ++i) {
+        for (int j = 0; j <= i; ++j) {
+          sum2[k++] = x[i] * x[j];
+        }
+      }
+
+    } else {
+      int k = 0;
+      for (int i = 0; i < dimension; ++i) {
+
+        if (x[i] < min[i]) {
+          min[i]        = x[i];
+          minIndices[i] = n;
+        } else if (x[i] > max[i]) {
+          max[i]        = x[i];
+          maxIndices[i] = n;
+        }
+
+        sum[i] += x[i];
+        for (int j = 0; j <= i; ++j) {
+          sum2[k++] += x[i] * x[j];
+        }
+
+      }
+    }
+
+    ++n;
+
+  }
+
+  /** Add all points of an array to the instance.
+   * @param points array of points
+   * @exception IllegalArgumentException if there is a dimension
+   * mismatch between these points and the ones already added (this
+   * cannot happen when the instance is empty)
+   */
+  public void add(double[][] points) {
+    for (int i = 0; i < points.length; ++i) {
+      add(points[i]);
+    }
+  }
+
+  /** Add all the points of another sample to the instance.
+   * @param s samples to add
+   * @exception IllegalArgumentException if there is a dimension
+   * mismatch between this sample points and the ones already added
+   * (this cannot happen when the instance is empty)
+   */
+  public void add(VectorialSampleStatistics s) {
+
+    if (s.n == 0) {
+      // nothing to add
+      return;
+    }
+
+    if (n == 0) {
+
+      dimension = s.dimension;
+      allocate();
+
+      System.arraycopy(s.min,         0, min,         0, dimension);
+      System.arraycopy(s.minIndices,  0, minIndices,  0, dimension);
+      System.arraycopy(s.max,         0, max,         0, dimension);
+      System.arraycopy(s.maxIndices,  0, maxIndices,  0, dimension);
+      System.arraycopy(s.sum,         0, sum,         0, dimension);
+      System.arraycopy(s.sum2,        0, sum2,        0, sum2.length);
+
+    } else {
+      int k = 0;
+
+      for (int i = 0; i < dimension; ++i) {
+
+        if (s.min[i] < min[i]) {
+          min[i]        = s.min[i];
+          minIndices[i] = n;
+        } else if (s.max[i] > max[i]) {
+          max[i]        = s.max[i];
+          maxIndices[i] = n;
+        }
+
+        sum[i] += s.sum[i];
+        for (int j = 0; j <= i; ++j) {
+          sum2[k] += s.sum2[k];
+          ++k;
+        }
+
+      }
+
+    }
+
+    n += s.n;
+
+  }
+
+  /** Get the number of points in the sample.
+   * @return number of points in the sample
+   */
+  public int size() {
+    return n;
+  }
+
+  /** Get the minimal value in the sample.
+   * <p>Since all components of the sample vector can reach their
+   * minimal value at different times, this vector should be
+   * considered as gathering all minimas of all components. The index
+   * of the sample at which the minimum was encountered can be
+   * retrieved with the {@link #getMinIndices getMinIndices}
+   * method.</p>
+   * @return minimal value in the sample (the array is a reference to
+   * an internal array that changes each time something is added to
+   * the instance, the caller should neither change it nor rely on its
+   * value in the long term)
+   * @see #getMinIndices
+   */
+  public double[] getMin() {
+    return min;
+  }
+
+  /** Get the indices at which the minimal value occurred in the sample.
+   * @return a vector reporting at which occurrence each component of
+   * the sample reached its minimal value (the array is a reference to
+   * an internal array that changes each time something is added to
+   * the instance, the caller should neither change it nor rely on its
+   * value in the long term)
+   * @see #getMin
+   */
+  public int[] getMinIndices() {
+    return minIndices;
+  }
+
+  /** Get the maximal value in the sample.
+   * <p>Since all components of the sample vector can reach their
+   * maximal value at different times, this vector should be
+   * considered as gathering all maximas of all components. The index
+   * of the sample at which the maximum was encountered can be
+   * retrieved with the {@link #getMaxIndices getMaxIndices}
+   * method.</p>
+   * @return maximal value in the sample (the array is a reference to
+   * an internal array that changes each time something is added to
+   * the instance, the caller should neither change it nor rely on its
+   * value in the long term)
+   * @see #getMaxIndices
+   */
+  public double[] getMax() {
+    return max;
+  }
+
+  /** Get the indices at which the maximal value occurred in the sample.
+   * @return a vector reporting at which occurrence each component of
+   * the sample reached its maximal value (the array is a reference to
+   * an internal array that changes each time something is added to
+   * the instance, the caller should neither change it nor rely on its
+   * value in the long term)
+   * @see #getMax
+   */
+  public int[] getMaxIndices() {
+    return maxIndices;
+  }
+
+  /** Get the mean value of the sample.
+   * @param mean placeholder where to store the array, if null a new
+   * array will be allocated
+   * @return mean value of the sample or null if the sample is empty
+   * and hence the dimension of the vectors is still unknown
+   * (reference to mean if it was non-null, reference to a new array
+   * otherwise)
+   */
+  public double[] getMean(double[] mean) {
+    if (n == 0) {
+      return null;
+    }
+    if (mean == null) {
+      mean = new double[dimension];
+    }
+    for (int i = 0; i < dimension; ++i) {
+      mean[i] = sum[i] / n;
+    }
+    return mean;
+  }
+
+  /** Get the covariance matrix of the underlying law.
+   * This method estimate the covariance matrix considering that the
+   * data available are only a <em>sample</em> of all possible
+   * values. This value is the sample covariance matrix (as opposed
+   * to the population covariance matrix).
+   * @param covariance placeholder where to store the matrix, if null
+   * a new matrix will be allocated
+   * @return covariance matrix of the underlying or null if the
+   * sample has less than 2 points (reference to covariance if it was
+   * non-null, reference to a new matrix otherwise)
+   */
+  public SymetricalMatrix getCovarianceMatrix(SymetricalMatrix covariance) {
+
+    if (n < 2) {
+      return null;
+    }
+
+    if (covariance == null) {
+      covariance = new SymetricalMatrix(dimension);
+    }
+
+    double c = 1.0 / (n * (n - 1));
+    int k = 0;
+    for (int i = 0; i < dimension; ++i) {
+      for (int j = 0; j <= i; ++j) {
+        double e = c * (n * sum2[k] - sum[i] * sum[j]);
+        covariance.setElementAndSymetricalElement(i, j, e);
+        ++k;
+      }
+    }
+
+    return covariance;
+
+    }
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/VectorialSampleStatistics.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/doc-files/org_spaceroots_mantissa_random_classes.png
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/doc-files/org_spaceroots_mantissa_random_classes.png?view=auto&rev=476930
==============================================================================
Binary file - no diff available.

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/doc-files/org_spaceroots_mantissa_random_classes.png
------------------------------------------------------------------------------
    svn:mime-type = image/png

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/package.html
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/package.html?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/package.html (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/package.html Sun Nov 19 13:36:15 2006
@@ -0,0 +1,33 @@
+<html>
+<body>
+This package provides classes to perform some random draws and
+statistical anaysis.
+
+<p>The aim of this package is to provide the basic components needed
+to generate random variables (correlated or not in the case of
+vectorial variables) that could be used in a simulation application
+and to provide some basic statistical classes in order to analyze the
+simulation results.</p>
+
+<p>Vectorial generators are build by embedding a {@link
+org.spaceroots.mantissa.random.NormalizedRandomGenerator normalized}
+scalar generator into the vectorial generator classes {@link
+org.spaceroots.mantissa.random.CorrelatedRandomVectorGenerator
+CorrelatedRandomVectorGenerator} or {@link
+org.spaceroots.mantissa.random.UncorrelatedRandomVectorGenerator
+UncorrelatedRandomVectorGenerator} that will be responsible for
+packaging all numbers into vectors with the specified mean values, standard
+deviations and correlation coefficients. Since most practical problems
+make the assumption the probability distribution is a gaussian one,
+the normalized generator will often be an instance of {@link
+org.spaceroots.mantissa.random.GaussianRandomGenerator
+GaussianRandomGenerator}, but uniform distribution are also available
+using instances of {@link
+org.spaceroots.mantissa.random.UniformRandomGenerator
+UniformRandomGenerator}.</p>
+
+<img src="doc-files/org_spaceroots_mantissa_random_classes.png" />
+
+@author L. Maisonobe
+</body>
+</html>

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/random/package.html
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/roots/BrentSolver.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/roots/BrentSolver.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/roots/BrentSolver.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/roots/BrentSolver.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,188 @@
+// 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.spaceroots.mantissa.roots;
+
+import org.spaceroots.mantissa.functions.scalar.ComputableFunction;
+import org.spaceroots.mantissa.functions.FunctionException;
+
+/** This class implements the Brent algorithm to compute the roots of
+ * a function in an interval.
+
+ * This class is basically a translation in Java of a fortran
+ * implementation found at netlib (<a
+ * href="http://www.netlib.org/fmm/zeroin.f">zeroin.f</a>).
+
+ * @version $Id: BrentSolver.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class BrentSolver implements RootsFinder {
+
+  /** IEEE 754 epsilon . */
+  private static final double epsilon = Math.pow(2.0, -52);
+
+  /** Root found. */
+  private double root;
+
+  /** Simple constructor.
+   * Build a Brent solver
+   */
+  public BrentSolver() {
+    root = Double.NaN;
+  }
+
+  /** Solve a function in a given interval known to contain a root.
+   * @param function function for which a root should be found
+   * @param checker checker for the convergence of the function
+   * @param maxIter maximal number of iteration allowed
+   * @param x0 abscissa of the lower bound of the interval
+   * @param f0 value of the function the lower bound of the interval
+   * @param x1 abscissa of the higher bound of the interval
+   * @param f1 value of the function the higher bound of the interval
+   * @return true if a root has been found in the given interval
+   */
+  public boolean findRoot(ComputableFunction function,
+                          ConvergenceChecker checker,
+                          int maxIter,
+                          double x0, double f0, double x1, double f1)
+    throws FunctionException {
+
+    double a  = x0;
+    double fa = f0;
+    double b  = x1;
+    double fb = f1;
+
+    double c  = a;
+    double fc = fa;
+
+    double d  = b - a;
+    double e  = d;
+
+    double tolS;
+    for (int iter = 0; iter < maxIter; ++iter) {
+
+      if (Math.abs(fc) < Math.abs(fb)) {
+        // invert points
+        a  = b;
+        b  = c;
+        c  = a;
+        fa = fb;
+        fb = fc;
+        fc = fa;
+      }
+
+      tolS = 2 * epsilon * Math.abs(b);
+      double xm = 0.5 * (c - b);
+
+      // convergence test
+      double xLow, fLow, xHigh, fHigh;
+      if (b < c) {
+        xLow   = b;
+        fLow   = fb;
+        xHigh  = c;
+        fHigh  = fc;
+      } else {
+        xLow   = c;
+        fLow   = fc;
+        xHigh  = b;
+        fHigh  = fb;
+      }
+
+      switch (checker.converged(xLow, fLow, xHigh, fHigh)) {
+      case ConvergenceChecker.LOW :
+        root = xLow;
+        return true;
+      case ConvergenceChecker.HIGH :
+        root = xHigh;
+        return true;
+      default :
+        if ((Math.abs(xm) < tolS) || (Math.abs(fb) < Double.MIN_VALUE)) {
+          root = b;
+          return true;
+        }
+      }
+
+      if ((Math.abs(e) < tolS) || (Math.abs(fa) <= Math.abs(fb))) {
+        // use bisection method
+        d = xm;
+        e = d;
+      } else {
+        // use secant method
+        double p, q, r, s;
+        s = fb / fa;
+        if (Math.abs(a - c) < epsilon * Math.max(Math.abs(a), Math.abs(c))) {
+          // linear interpolation using only b and c points
+          p = 2.0 * xm * s;
+          q = 1.0 - s;
+        } else {
+          // inverse quadratic interpolation using a, b and c points
+          q = fa / fc;
+          r = fb / fc;
+          p = s * (2.0 * xm * q * (q - r) - (b - a) * (r - 1.0));
+          q = (q - 1.0) * (r - 1.0) * (s - 1.0);
+        }
+
+        // signs adjustment
+        if (p > 0.0) {
+          q = -q;
+        } else {
+          p = -p;
+        }
+
+        // is interpolation acceptable ?
+        if (((2.0 * p) < (3.0 * xm * q - Math.abs(tolS * q)))
+            &&
+            (p < Math.abs(0.5 * e * q))) {
+          e = d;
+          d = p / q;
+        } else {
+          // no, we need to fall back to bisection
+          d = xm;
+          e = d;
+        }
+      }
+
+      // complete step
+      a  = b;
+      fa = fb;
+      b += ((Math.abs(d) > tolS) ? d : (xm > 0.0 ? tolS : -tolS));
+      fb = function.valueAt(b);
+
+      if (fb * fc > 0) {
+        c  = a;
+        fc = fa;
+        d  = b - a;
+        e  = d;
+      }
+
+    }
+
+    // we have exceeded the maximal number of iterations
+    return false;
+
+  }
+
+  /** Get the abscissa of the root.
+   * @return abscissa of the root
+   */
+  public double getRoot() {
+    return root;
+  }
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/roots/BrentSolver.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/roots/ConvergenceChecker.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/roots/ConvergenceChecker.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/roots/ConvergenceChecker.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/roots/ConvergenceChecker.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,57 @@
+// 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.spaceroots.mantissa.roots;
+
+/** This interface specifies methods to check if a root-finding
+ * algorithm has converged.
+
+ * Deciding if convergence has been reached is a problem-dependent
+ * issue. The user should provide a class implementing this interface
+ * to allow the root-finding algorithm to stop its search according to
+ * the problem at hand.
+
+ * @version $Id: ConvergenceChecker.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public interface ConvergenceChecker {
+
+  /** Indicator for no convergence. */
+  public static final int NONE = 0;
+
+  /** Indicator for convergence on the lower bound of the interval. */
+  public static final int LOW  = 1;
+
+  /** Indicator for convergence on the higher bound of the interval. */
+  public static final int HIGH = 2;
+
+  /** Check if the root-finding algorithm has converged on the interval.
+   * The interval defined by the arguments contains one root (if there
+   * was at least one in the initial interval given by the user to the
+   * root-finding algorithm, of course)
+   * @param xLow abscissa of the lower bound of the interval
+   * @param fLow value of the function the lower bound of the interval
+   * @param xHigh abscissa of the higher bound of the interval
+   * @param fHigh value of the function the higher bound of the interval
+   * @return convergence indicator, must be one of {@link #NONE},
+   * {@link #LOW} or {@link #HIGH}
+   */
+  public int converged (double xLow, double fLow, double xHigh, double fHigh);
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/roots/ConvergenceChecker.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/roots/RootsFinder.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/roots/RootsFinder.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/roots/RootsFinder.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/roots/RootsFinder.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,54 @@
+// 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.spaceroots.mantissa.roots;
+
+import org.spaceroots.mantissa.functions.scalar.ComputableFunction;
+import org.spaceroots.mantissa.functions.FunctionException;
+
+/** This interface specifies root-finding methods for scalar
+ * functions.
+
+ * @version $Id: RootsFinder.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public interface RootsFinder {
+
+  /** Solve a function in a given interval known to contain a root.
+   * @param function function for which a root should be found
+   * @param checker checker for the convergence of the function
+   * @param maxIter maximal number of iteration allowed
+   * @param x0 abscissa of the lower bound of the interval
+   * @param f0 value of the function the lower bound of the interval
+   * @param x1 abscissa of the higher bound of the interval
+   * @param f1 value of the function the higher bound of the interval
+   * @return true if a root has been found in the given interval
+   */
+  public boolean findRoot(ComputableFunction function,
+                          ConvergenceChecker checker,
+                          int maxIter,
+                          double x0, double f0, double x1, double f1)
+    throws FunctionException;
+
+  /** Get the abscissa of the root.
+   * @return abscissa of the root
+   */
+  public double getRoot();
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/roots/RootsFinder.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/utilities/ArrayMapper.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/utilities/ArrayMapper.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/utilities/ArrayMapper.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/utilities/ArrayMapper.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,141 @@
+// 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.spaceroots.mantissa.utilities;
+
+import java.util.ArrayList;
+import java.util.Iterator;
+
+/**
+ * This class dispatch data between an array and several domain objects.
+
+ * This class handles all the burden of mapping each domain object it
+ * handles to a slice of a single array.
+
+ * @see ArraySliceMappable
+
+ * @version $Id: ArrayMapper.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class ArrayMapper {
+
+  /** Simple constructor.
+   * Build an empty array mapper
+   */
+  public ArrayMapper() {
+    domainObjects = new ArrayList();
+    size          = 0;
+    internalData  = null;
+  }
+
+  /** Simple constructor.
+   * Build an array mapper managing one object. Other objects can be
+   * added later using the {@link #manageMappable manageMappable}
+   * method. This call is equivalent to build the mapper with the
+   * default constructor and adding the object.
+   * @param object domain object to handle
+   */
+  public ArrayMapper(ArraySliceMappable object) {
+
+    domainObjects = new ArrayList();
+    domainObjects.add(new ArrayMapperEntry(object, 0));
+
+    size = object.getStateDimension();
+
+    internalData = new double [size];
+
+  }
+
+  /** Take a new domain object into account.
+   * @param object domain object to handle
+   */
+  public void manageMappable(ArraySliceMappable object) {
+
+    domainObjects.add(new ArrayMapperEntry(object, size));
+
+    size += object.getStateDimension();
+
+    if (internalData != null) {
+      internalData = new double [size];
+    }
+
+  }
+
+  /** Get the internal data array.
+   * @return internal data array
+   */
+  public double[] getInternalDataArray() {
+    if (internalData == null) {
+      internalData = new double [size];
+    }
+    return internalData;
+  }
+
+  /** Map data from the internal array to the domain objects.
+   */
+  public void updateObjects() {
+    if (internalData == null) {
+      internalData = new double [size];
+    }
+    updateObjects(internalData);
+  }
+
+  /** Map data from the specified array to the domain objects.
+   * @param data flat array holding the data to dispatch
+   */
+  public void updateObjects(double[] data) {
+    for (Iterator iter = domainObjects.iterator(); iter.hasNext();) {
+      ArrayMapperEntry entry = (ArrayMapperEntry) iter.next();
+      entry.object.mapStateFromArray(entry.offset, data);
+    }
+  }
+
+  /** Map data from the domain objects to the internal array.
+   */
+  public void updateArray() {
+    if (internalData == null) {
+      internalData = new double [size];
+    }
+    updateArray(internalData);
+  }
+
+  /** Map data from the domain objects to the specified array.
+   * @param data flat array where to put the data
+   */
+  public void updateArray(double[] data) {
+    for (Iterator iter = domainObjects.iterator(); iter.hasNext();) {
+      ArrayMapperEntry entry = (ArrayMapperEntry) iter.next();
+      entry.object.mapStateToArray(entry.offset, data);
+    }
+  }
+
+  /** Container for all handled objects. */
+  private ArrayList domainObjects;
+
+  /** Total number of scalar elements handled.
+   * (size of the array)
+   */
+  private int size;
+
+  /** Flat array holding all data.
+   * This is null as long as nobody uses it (lazy creation)
+   */
+  private double[] internalData;
+
+}

Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/utilities/ArrayMapper.java
------------------------------------------------------------------------------
    svn:eol-style = native



---------------------------------------------------------------------
To unsubscribe, e-mail: commons-dev-unsubscribe@jakarta.apache.org
For additional commands, e-mail: commons-dev-help@jakarta.apache.org


Mime
View raw message