gov.sandia.cognition.math.LentzMethod Maven / Gradle / Ivy
/*
* File: LentzMethod.java
* Authors: Kevin R. Dixon
* Company: Sandia National Laboratories
* Project: Cognitive Foundry
*
* Copyright Apr 16, 2008, Sandia Corporation.
* Under the terms of Contract DE-AC04-94AL85000, there is a non-exclusive
* license for use of this work by or on behalf of the U.S. Government.
* Export of this program may require a license from the United States
* Government. See CopyrightHistory.txt for complete details.
*
*/
package gov.sandia.cognition.math;
import gov.sandia.cognition.algorithm.AbstractAnytimeAlgorithm;
import gov.sandia.cognition.annotation.PublicationReference;
import gov.sandia.cognition.annotation.PublicationType;
/**
* This class implements Lentz's method for evaluating continued fractions.
* This is an iterative, approximate approach that has good convergence. By
* specifying the fraction as f(x) = b0 + a1/(b1+a1/(b2+a3/(b3+...))).
* Each iteration, we take the (a,b) coefficients and improve the approximation
* of f(x) and terminate after too many iterations or convergence.
*
* @author Kevin R. Dixon
* @since 2.1
*/
@PublicationReference(
author={
"William H. Press",
"Saul A. Teukolsky",
"William T. Vetterling",
"Brian P. Flannery"
},
title="Numerical Recipes in C, Second Edition",
type=PublicationType.Book,
year=1992,
pages={169,171},
url="http://www.nrbook.com/a/bookcpdf.php"
)
public class LentzMethod
extends AbstractAnytimeAlgorithm
{
/**
* Default value to keep denominators from equaling 0.0, default {@value}
*/
final static double DEFAULT_MIN_DENOMINATOR = 1e-30;
/**
* Default max iterations, {@value}
*/
final static int DEFAULT_MAX_ITERATIONS = 100;
/**
* Default convergence tolerance, {@value}
*/
final static double DEFAULT_TOLERANCE = 3e-7;
/**
* Tolerance of the algorithm for convergence
*/
private double tolerance;
/**
* Value to keep denominators from equaling 0.0
*/
private double minDenominator;
/**
* Current value of the "C" variable in Lentz's method
*/
private double currentC;
/**
* Current value of the "D" variable in Lentz's method
*/
private double currentD;
/**
* Value of the continued fraction, null if not valid
*/
private Double result;
/**
* Running (intermediate) value of the fraction value
*/
private double fractionValue;
/**
* Flag to keep going or stop
*/
private boolean keepGoing;
/**
* Creates a new instance of LentzMethod
*/
public LentzMethod()
{
this( DEFAULT_MAX_ITERATIONS, DEFAULT_TOLERANCE, DEFAULT_MIN_DENOMINATOR );
}
/**
* Creates a new instance of LentzMethod
*
* @param maxIterations
* Maximum number of iterations before stopping
* @param tolerance
* Tolerance of the algorithm for convergence
* @param minDenominator
* Value to keep denominators from equaling 0.0
*/
public LentzMethod(
int maxIterations,
double tolerance,
double minDenominator )
{
super( maxIterations );
this.setTolerance(tolerance);
this.setMinDenominator(minDenominator);
this.keepGoing = false;
this.result = null;
}
/**
* Initializes Lentz's method using the given values
* @param b0
* Initial "b" coefficient of the continued fraction
* @return
* keep going flag
*/
public boolean initializeAlgorithm(
double b0 )
{
this.fireAlgorithmStarted();
result = null;
this.fractionValue = b0;
if( Math.abs(this.fractionValue) < this.minDenominator )
{
this.fractionValue = this.minDenominator;
}
this.currentC = this.fractionValue;
this.currentD = 0.0;
this.iteration = 0;
this.keepGoing = (this.iteration < this.maxIterations);
return this.keepGoing;
}
/**
* Step of Lentz's iteration
* @param a
* coefficient "a" from the next recursion of the continued fraction
* @param b
* coefficient "b" from the next recursion of the continued fraction
* @return
* keep going flag
*/
public boolean iterate(
double a,
double b )
{
// If we should have already stopped, but you want us to keep going,
// then throw an exception
if( this.keepGoing == false )
{
throw new OperationNotConvergedException(
"Trying to iterate when keepGoing is false!" );
}
this.fireStepStarted();
this.iteration++;
this.currentD = b + a * this.currentD;
this.currentC = b + a / this.currentC;
if( Math.abs(this.currentD) < this.minDenominator )
{
this.currentD = this.minDenominator;
}
this.currentD = 1.0 / this.currentD;
if( Math.abs(this.currentC) < this.minDenominator )
{
this.currentC = this.minDenominator;
}
double delta = this.currentC * this.currentD;
this.fractionValue *= delta;
// We're done!
if( Math.abs(delta-1.0) < this.tolerance )
{
this.result = this.getFractionValue();
this.keepGoing = false;
}
// We've iterated too much, so stop
if( this.iteration >= this.maxIterations )
{
this.keepGoing = false;
}
this.fireStepEnded();
// Fire the algorithm-ended message, if appropriate
if( this.keepGoing == false )
{
this.fireAlgorithmEnded();
}
return this.keepGoing;
}
/**
* Getter for tolerance
* @return
* Tolerance of the algorithm for convergence
*/
public double getTolerance()
{
return this.tolerance;
}
/**
* Setter for tolerance
* @param tolerance
* Tolerance of the algorithm for convergence
*/
public void setTolerance(
double tolerance )
{
this.tolerance = tolerance;
}
/**
* Getter for minDenominator
*
* @return The minimum denominator.
*/
public double getMinDenominator()
{
return this.minDenominator;
}
/**
* Setter for minDenominator
* @param minDenominator
* Value to keep denominators from equaling 0.0
*/
public void setMinDenominator(
double minDenominator )
{
this.minDenominator = minDenominator;
}
public Double getResult()
{
return this.result;
}
public void stop()
{
this.keepGoing = false;
}
/**
* Getter for keepGoing
* @return
* Flag to keep going or stop
*/
public boolean getKeepGoing()
{
return this.keepGoing;
}
/**
* Getter for fractionValue
* @return
* Running (intermediate) value of the fraction value
*/
public double getFractionValue()
{
return this.fractionValue;
}
/**
* Setter for fractionValue
* @param fractionValue
* Running (intermediate) value of the fraction value
*/
protected void setFractionValue(
double fractionValue )
{
this.fractionValue = fractionValue;
}
}