Implmentation of contnued fraction

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

Hi,

I’m not sure, about your question.

  • 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