The big_gamma
function,
,
is given in
terms of a continued fraction expansion:
= | |||
= | (2.27) |
Unlike a power series expansion the continued fraction expansion seemingly cannot be evaluated in steps from ``left to right'', i.e., from the lowest terms, gradually adding higher and higher terms to the sum. Here it looks like we have to start from the end, sic! That is, we have to guess the last term, and then gradually build up the whole fraction starting from the last term in the expansion and moving towards the beginning.
But it can be done the right way, i.e., from ``left to right''.
Consider the following generic continued fraction expansions:
f_{n-1}(x) | = | ||
= | |||
f_{n}(x) | = | ||
= | (2.28) |
We can help ourselves with Maxima in this task.
First we
define f_{1} and extract its numerator and denominator:
biga[2]
in a way that shows the presence of biga[1]
in it!
Let us now evaluate f_{3}:
ratsubst
, which takes three arguments: the name of a new
variable, its value, and the expression into which the variable
should be substituted:
A_{3} = b_{3} A_{2} + a_{3} A_{1} | (2.29) |
B_{3} = b_{3} B_{2} + a_{3} B_{1} | (2.30) |
The evaluation of a continued fraction expansion now looks much easier. All that we have to do is to
Looking once again at the continued fraction expansion for
:
Observe that a regular sequencing begins only with a_{2} and b_{2}, not with a_{1} and b_{1}. In fact the first term in the continued expansion for is the odd one. Consequently we shall begin our iterations from evaluating A_{2} and B_{2}, and we'll use A_{0}, B_{0}, A_{1}, and B_{1} as our startup values.
Using mathematicised notation our algorithm is going to be as follows:
assume | |
DO | |
now a_{n} becomes a_{2} | |
now b_{n} becomes b_{2} | |
now A_{0} becomes A_{2} | |
now B_{0} becomes B_{2} | |
now a_{n} becomes a_{3} | |
now b_{n} becomes b_{3} | |
now A_{1} becomes A_{3} | |
now B_{1} becomes B_{3} | |
UNTIL |
At this stage we have carefully mapped the algorithm, without
attempting any economies. But having a closer look at the DO
loop as it's been outlined above shows that we can make some savings.
First, we don't really need to maintain the variables a_{n} and b_{n}. Instead
we can use n-a, 1, n and x directly:
DO | |
now A_{0} becomes A_{2} | |
now B_{0} becomes B_{2} | |
now A_{1} becomes A_{3} | |
now B_{1} becomes B_{3} | |
This saves us 4 assignments within the loop. Assignments are not as cheap as they seem to be, because data movement is involved
Figure 2.13 shows subprogram big_gamma
. You should
be able to identify clearly all elements. The algorithm is exactly
as specified above.
Having now implemented all required functions we can not only calculate the
goodness of fit in our problem, which happens to be Q=0.053, but we can
test the implementation of function goodness
itself.
Figure 2.14 shows a short program, which calculates 1 - Q(a, x) = P(a, x), and writes the results on a data file, try.out. This file can then be displayed using Gnuplot and the graph compared against similar graphs found in the literature.
Figure 2.15 shows the Gnuplot command file used to convert try.out into a graph shown in Figure 2.16.