# Implmentation of contnued fraction

Hi all,
How can we implement continued fraction statements in fortran?
For example:

Hi,

• Do you want, the value of “f” from the two lists of integers a_i and b_i?
=> it is relatively easy.
• or do you want the decomposition of “f” into the two lists of integers?
=> You need an algorithm for that (see wikipedia)

By the way, for continued fraction, normally, the a_i are equal to one.

1 Like

The question is:

As stated by @gardhor, it is relatively easy if you use the Euler-Wallis recurrence relations. Then you just need to write a recursive function:

module continued_fractions

implicit none
integer, parameter :: pi_seq(9) = [3, 7, 15, 1, 292, 1, 1, 1, 2]
!! Simple continued fraction expansion of Pi taken from https://oeis.org/A001203

contains

!>  Returns the *simple* continued fraction of the sequence b
recursive function cf(b) result(pq)
integer, intent(in) :: b(:)
real :: pq(2)

real :: pq_n1(2), pq_n2(2)
integer :: n

n = size(b)
if (n > 2) then
pq_n2 = cf(b(1:n-2))
pq_n1 = cf(b(1:n-1))
pq = b(n)*pq_n1 + pq_n2
else if (n == 2) then
pq_n1 = cf(b(1:1))
pq(1) = b(2)*pq_n1(1) + 1.0
pq(2) = b(2)
else if (n == 1) then
pq(1) = b(1)
pq(2) = 1.0
end if
end function

end module

program main

use continued_fractions, only: cf, pi_seq
implicit none

real :: f, pq(2)
real, parameter :: pi = 4.0*atan(1.0)

pq = cf(pi_seq)
f = pq(1)/pq(2)

print *, "continued fraction:", f
print *, "true value:", pi

end program


I did not implement the scaling part to prevent overflow. The reason you need it however, is the values of p_n and q_n quickly become large integers as you add more terms. In the example above with just 9 terms of the continued fractions, the values in pq are:

833719.000       265381.000


To reach the full 20 fractions, it seems like you will need the rescaling.

1 Like

As a cool side note, the fourth convergent of \pi, [3;7,15,1] = 355/113, is also known as Milü after a Chinese astronomer born 429 AD. It is within 0.000009% of the value of \pi.

5 Likes