2
\$\begingroup\$

I have the following recursive function to compute the continued fraction:

$$s_i = \frac{a_{i}}{b_{i} + s_{i+1}}$$

The relative error is defined as:

$$\text{Delta} = \frac{s_n - s_{n-1}}{s_{n}}$$

Here is the implementation in Python, of Lentz's method for computing continued fractions:

import numpy as np

def n_continued_fraction(a, b, err, i_min=3):
    f = lambda i : a(i)
    g = lambda i : b(i)
    d = lambda i : 1./g(0) if i == 0 else 1./(g(i) + f(i) * d(i-1))
    h = lambda i : 0 if i == -1 else f(0)/g(0) if i == 0 else (h(i-1) * g(i) + h(i-2) * d(i-1) * f(i)) * d(i)
    delta_h = lambda i : h(i) - h(i-1)

    i_final = i_min

    while np.abs(delta_h(i_final)) > err * np.abs(h(i_final)):
        i_final += 1

    return h(i_final)

and you can test it with e.g.

a = lambda i : (i+1.)
b = lambda i : (i+5.)

n_continued_fraction(a, b, 1.e-16)

Is there a better method to do it? Am I doing anything wrong?

\$\endgroup\$
2
\$\begingroup\$

Nice code, it does what it should and it is short. But it can be improved :)


General python tips:

  • Use type annotation. You can look at the type annotation of Sam Stafford in the previous answer, it will make your code more readable and will prevent bugs.
  • Remove redundant imports. Is numpy really needed here? You can use a simple abs.
  • Don't use lambda functions. According to PEP 8 "do not assign lambda expression, use a def", there are several reasons why, mostly due to readability.

Optimisation tips:
Note these tips contradict each other. Choose one and use it but you cannot use both.

  • If N is the number of elements until convergence, and a and b is the calculation time of each function (Assuming equal for every x) The time complexity now is O(N ^ 2 * a * b), for each i you calculate the entire series over and over.
    If you will use exponential search you can reduce the complexity to O(N * log N * a * b). Which will be much much faster!
  • Readability vs. Runtime. If the readability is important to you you can ignore this tip. You can replace the recursion with a dynamic programing approach, this will remove the recursion and the iteration over i, the time complexity will be O(N * a * b) but you will need to rewrite a lot of the code.
\$\endgroup\$
  • \$\begingroup\$ Thanks a lot. Can I get a hand with the second tip, i.e., replace the recursion with dynamic programming? \$\endgroup\$ – John Jan 12 at 15:16
  • \$\begingroup\$ Well I would like to help you but you will need to answer Peilonrayz question in the comments first. @user216427 \$\endgroup\$ – Yonlif Jan 13 at 18:54
  • \$\begingroup\$ Sure, I've tried to answer them as best as I could. \$\endgroup\$ – John Jan 13 at 22:36
  • \$\begingroup\$ This function will be used a few hundred thousand times, so I am interested in achieving the fastest performance possible, even if I have to sacrifice readability. \$\endgroup\$ – John Jan 14 at 13:22
  • \$\begingroup\$ I'm interested in how you'd get exponential search to work with recursive sequence with \$O(n\log{n})\$ time. I don't know why you're saying DP is hard to read. It's fairly simple, and as readable as the OP's. \$\endgroup\$ – Peilonrayz Jan 14 at 16:39
0
\$\begingroup\$

Since lambda expressions are hard to add types and decorators to, I'd start off by writing it like this:

import numpy as np
from typing import Callable

def n_continued_fraction(
    a: Callable[[float], float], 
    b: Callable[[float], float], 
    err: float, 
    i_min: float = 3.0
) -> float:

    def d(i: float) -> float:
        if i == 0:
            return 1.0 / b(0)
        return 1.0 / (b(i) + a(i) * d(i - 1))

    def h(i: float) -> float:
        if i == -1:
            return 0
        if i == 0:
            return a(0) / b(0)
        return (h(i - 1) * b(i) + h(i - 2) * d(i - 1) * a(i)) * d(i)

    def delta_h(i: float) -> float:
        return h(i) - h(i - 1)

    i = i_min
    while np.abs(delta_h(i)) > err * np.abs(h(i)):
        i += 1
    return h(i)

In terms of performance, note that Python doesn't do tail recursion optimization, so if the stack gets deep enough the script will exit.

However! Defining functions with def opens up the possibility of using the power of decorators to memoize those inner functions, which might help you get around the recursion limitations; I'll leave experimenting with that as an exercise for the reader. :) https://www.geeksforgeeks.org/memoization-using-decorators-in-python/

\$\endgroup\$
  • \$\begingroup\$ How's this better than the OPs code? \$\endgroup\$ – Peilonrayz Jan 10 at 22:32
  • \$\begingroup\$ i.imgur.com/bJQPWcy.png \$\endgroup\$ – Sam Stafford Jan 10 at 22:34
  • 1
    \$\begingroup\$ You're so funny. I can't contain my laughter. \$\endgroup\$ – Peilonrayz Jan 10 at 22:35
  • \$\begingroup\$ How "deep" does the stack need to get before it exits? Is this completely resolved using separate functions? \$\endgroup\$ – John Jan 10 at 22:41
  • \$\begingroup\$ It's somewhere on the order of thousands. Using separate functions doesn't help if one calls the other before returning; you can think of stack depth as "how many function calls are still waiting to return", total across the program. Memoization can potentially help by letting some of those calls return earlier than they would have otherwise. \$\endgroup\$ – Sam Stafford Jan 11 at 0:12

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.