KickJava   Java API By Example, From Geeks To Geeks.

Java > Open Source Codes > org > apache > commons > math > util > ContinuedFraction


1 /*
2  * Copyright 2003-2004 The Apache Software Foundation.
3  *
4  * Licensed under the Apache License, Version 2.0 (the "License");
5  * you may not use this file except in compliance with the License.
6  * You may obtain a copy of the License at
7  *
8  * http://www.apache.org/licenses/LICENSE-2.0
9  *
10  * Unless required by applicable law or agreed to in writing, software
11  * distributed under the License is distributed on an "AS IS" BASIS,
12  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13  * See the License for the specific language governing permissions and
14  * limitations under the License.
15  */
16 package org.apache.commons.math.util;
17
18 import java.io.Serializable JavaDoc;
19
20 import org.apache.commons.math.ConvergenceException;
21 import org.apache.commons.math.MathException;
22
23 /**
24  * Provides a generic means to evaluate continued fractions. Subclasses simply
25  * provided the a and b coefficients to evaluate the continued fraction.
26  *
27  * <p>
28  * References:
29  * <ul>
30  * <li><a HREF="http://mathworld.wolfram.com/ContinuedFraction.html">
31  * Continued Fraction</a></li>
32  * </ul>
33  * </p>
34  *
35  * @version $Revision$ $Date: 2005-02-26 05:11:52 -0800 (Sat, 26 Feb 2005) $
36  */
37 public abstract class ContinuedFraction implements Serializable JavaDoc {
38     
39     /** Serialization UID */
40     static final long serialVersionUID = 1768555336266158242L;
41     
42     /** Maximum allowed numerical error. */
43     private static final double DEFAULT_EPSILON = 10e-9;
44
45     /**
46      * Default constructor.
47      */
48     protected ContinuedFraction() {
49         super();
50     }
51
52     /**
53      * Access the n-th a coefficient of the continued fraction. Since a can be
54      * a function of the evaluation point, x, that is passed in as well.
55      * @param n the coefficient index to retrieve.
56      * @param x the evaluation point.
57      * @return the n-th a coefficient.
58      */
59     protected abstract double getA(int n, double x);
60
61     /**
62      * Access the n-th b coefficient of the continued fraction. Since b can be
63      * a function of the evaluation point, x, that is passed in as well.
64      * @param n the coefficient index to retrieve.
65      * @param x the evaluation point.
66      * @return the n-th b coefficient.
67      */
68     protected abstract double getB(int n, double x);
69
70     /**
71      * Evaluates the continued fraction at the value x.
72      * @param x the evaluation point.
73      * @return the value of the continued fraction evaluated at x.
74      * @throws MathException if the algorithm fails to converge.
75      */
76     public double evaluate(double x) throws MathException {
77         return evaluate(x, DEFAULT_EPSILON, Integer.MAX_VALUE);
78     }
79
80     /**
81      * Evaluates the continued fraction at the value x.
82      * @param x the evaluation point.
83      * @param epsilon maximum error allowed.
84      * @return the value of the continued fraction evaluated at x.
85      * @throws MathException if the algorithm fails to converge.
86      */
87     public double evaluate(double x, double epsilon) throws MathException {
88         return evaluate(x, epsilon, Integer.MAX_VALUE);
89     }
90
91     /**
92      * Evaluates the continued fraction at the value x.
93      * @param x the evaluation point.
94      * @param maxIterations maximum number of convergents
95      * @return the value of the continued fraction evaluated at x.
96      * @throws MathException if the algorithm fails to converge.
97      */
98     public double evaluate(double x, int maxIterations) throws MathException {
99         return evaluate(x, DEFAULT_EPSILON, maxIterations);
100     }
101
102     /**
103      * Evaluates the continued fraction at the value x.
104      *
105      * The implementation of this method is based on:
106      * <ul>
107      * <li>O. E-gecio-glu, C . K. Koc, J. Rifa i Coma,
108      * <a HREF="http://citeseer.ist.psu.edu/egecioglu91fast.html">
109      * On Fast Computation of Continued Fractions</a>, Computers Math. Applic.,
110      * 21(2--3), 1991, 167--169.</li>
111      * </ul>
112      *
113      * @param x the evaluation point.
114      * @param epsilon maximum error allowed.
115      * @param maxIterations maximum number of convergents
116      * @return the value of the continued fraction evaluated at x.
117      * @throws MathException if the algorithm fails to converge.
118      */
119     public double evaluate(double x, double epsilon, int maxIterations)
120         throws MathException
121     {
122         double[][] f = new double[2][2];
123         double[][] a = new double[2][2];
124         double[][] an = new double[2][2];
125
126         a[0][0] = getA(0, x);
127         a[0][1] = 1.0;
128         a[1][0] = 1.0;
129         a[1][1] = 0.0;
130
131         return evaluate(1, x, a, an, f, epsilon, maxIterations);
132     }
133
134     /**
135      * Evaluates the n-th convergent, fn = pn / qn, for this continued fraction
136      * at the value x.
137      * @param n the convergent to compute.
138      * @param x the evaluation point.
139      * @param a (n-1)-th convergent matrix. (Input)
140      * @param an the n-th coefficient matrix. (Output)
141      * @param f the n-th convergent matrix. (Output)
142      * @param epsilon maximum error allowed.
143      * @param maxIterations maximum number of convergents
144      * @return the value of the the n-th convergent for this continued fraction
145      * evaluated at x.
146      * @throws MathException if the algorithm fails to converge.
147      */
148     private double evaluate(
149         int n,
150         double x,
151         double[][] a,
152         double[][] an,
153         double[][] f,
154         double epsilon,
155         int maxIterations) throws MathException
156     {
157         double ret;
158
159         // create next matrix
160 an[0][0] = getA(n, x);
161         an[0][1] = 1.0;
162         an[1][0] = getB(n, x);
163         an[1][1] = 0.0;
164
165         // multiply a and an, save as f
166 f[0][0] = (a[0][0] * an[0][0]) + (a[0][1] * an[1][0]);
167         f[0][1] = (a[0][0] * an[0][1]) + (a[0][1] * an[1][1]);
168         f[1][0] = (a[1][0] * an[0][0]) + (a[1][1] * an[1][0]);
169         f[1][1] = (a[1][0] * an[0][1]) + (a[1][1] * an[1][1]);
170
171         // determine if we're close enough
172 if (Math.abs((f[0][0] * f[1][1]) - (f[1][0] * f[0][1])) <
173             Math.abs(epsilon * f[1][0] * f[1][1]))
174         {
175             ret = f[0][0] / f[1][0];
176         } else {
177             if (n >= maxIterations) {
178                 throw new ConvergenceException(
179                     "Continued fraction convergents failed to converge.");
180             }
181             // compute next
182 ret = evaluate(n + 1, x, f /* new a */
183             , an /* reuse an */
184             , a /* new f */
185             , epsilon, maxIterations);
186         }
187
188         return ret;
189     }
190 }
191
Popular Tags