Details
-
Improvement
-
Status: Closed
-
Major
-
Resolution: Fixed
-
2.0
-
None
-
None
Description
One can use a one-dimensional array (instead of Romberg's tableau) for extrapolating subsequent values.
Please have a look at following code fragments (which I've taken form the class RombergExtrapolator of
my MathLibrary). Feel free to use this code.
/**
- Default number of maximal extrapolation steps.
*/
public static int DEF_MAXIMAL_EXTRAPOLATION_COUNT = 8;
/**
- The approximation order. <br>
- Assume that f(h) is approximated by a function a(h), so that f(h) = a(h) +
- O(h<sup>p</sup>). We say that p is the approximation order.
*/
private int approximationOrder;
private int extrapolationCount = 0;
private double prevResult;
/**
- The estimate and tolerance may be used to deside wether to finalize the
- iteration process (|estimate| < tolerance).
*/
/** Holds the current estimated error. */
private double estimate;
/** Holds the current reached tolerance. */
private double tolerance;
private double result[] = new double[DEF_MAXIMAL_EXTRAPOLATION_COUNT + 1];;
/**
- Set the maximal number of subsequent extrapolation steps.
- @param maximalExtrapolationCount
- maximal extrapolation steps
*/
public void setMaximalExtrapolationCount(int maximalExtrapolationCount) { result = new double[maximalExtrapolationCount + 1]; }
/**
- Extrapolate a sequence of values by means of Romberg's algorithm.
- Therefore a polynomial of degree maximalExtraploationCount
- is used. Calculates the current estimate and tolerance using the
- approximation order.
- @param value
- value to extrapolate
- @return extrapolated value
*/
public double extrapolate(double value)
{
if (extrapolationCount == 0) { // first estimate estimate = value; tolerance = -1.0; prevResult = 0; }
int i, m, m1 = idx(extrapolationCount);
long k = (1 << approximationOrder);
int imin = Math.max(0, extrapolationCount - (result.length - 1));
result[m1] = value;
for (i = extrapolationCount - 1; i >= imin; i--)
{ m = idx(i); m1 = idx(i + 1); result[m] = (k * result[m1] - result[m]) / (k - 1); k <<= approximationOrder; } m1 = idx(i + 1);
estimate = result[m1] - prevResult;
tolerance = Math.abs(result[m1]) * relativeAccuracy + absoluteAccuracy;
prevResult = result[m1];
extrapolationCount++;
return result[m1];
}
/**
- Ring buffer index
*/
private int idx(int i) { return (i % result.length); }
