# What are some good ways to compute terms in recurrence relations?

**URL:** <https://discourse.julialang.org/t/what-are-some-good-ways-to-compute-terms-in-recurrence-relations/42016>\
**Category:** General Usage\
**Created:** [June 25, 2020, 12:18am UTC](https://discourse.julialang.org/t/what-are-some-good-ways-to-compute-terms-in-recurrence-relations/42016 "2020-06-25T00:18:40Z")\
**Posts on this page:** 10\
**Page:** 1

<div class="post-metadata">

**Author:** ![xiaodai](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/xiaodai/32/15937_2.png) [@xiaodai](https://discourse.julialang.org/u/xiaodai)\
**Post date:** [June 25, 2020, 12:18am UTC](https://discourse.julialang.org/t/what-are-some-good-ways-to-compute-terms-in-recurrence-relations/42016/1 "2020-06-25T00:18:40Z")

</div>

Doing a google search on “recurrence relations Julia” shows up MathToolkit.jl which has not been updated for a while and it’s got a neat function for detecting recurrence relations but not for solving them given the coefficients.

For recurrence relations that “reduces”, Memoize.jl works really well. E.g. Fibonaci sequence is just

```julia
using Memoize
@memoize function fib(n)
    if n <= 2
        return 1
    else
       return fib(n-1)+fib(n-2)
    end
end

```

And that’s it! However, I have issues with recurrence relations that don’t “reduce” e.g. if `f(n)` depends on `f(m)` where `m > n`.

E.g. say I am looking for the expected number of steps needed to reach finished line and the steps taken are controlled by a die where rolling 1 means take one step back and rolling 6 means go forward one step and rolling anything else means staying put. Then I get a recurrence relations like this

```julia
f(n) = 1 + f(n-1)//6 + f(n+1)//6 + f(n)*4//6

```

where `n` is the number of steps still to go in order to reach finish line.

I am just wondering if there’s a good way to solve this in Julia simply. I wrote some recursive function to solve it.

---

<div class="post-metadata">

**Author:** ![ericphanson](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ericphanson/32/215186_2.png) [@ericphanson](https://discourse.julialang.org/u/ericphanson)\
**Post date:** [June 25, 2020, 12:45am UTC](https://discourse.julialang.org/t/what-are-some-good-ways-to-compute-terms-in-recurrence-relations/42016/2 "2020-06-25T00:45:21Z")

</div>

This isn’t a general answer, but it seems like in the example there you could just solve for `f(n+1)` and then have a recurrence relation only in terms of lower values, at which point you can just memoize.

---

<div class="post-metadata">

**Author:** ![xiaodai](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/xiaodai/32/15937_2.png) [@xiaodai](https://discourse.julialang.org/u/xiaodai)\
**Post date:** [June 25, 2020, 12:58am UTC](https://discourse.julialang.org/t/what-are-some-good-ways-to-compute-terms-in-recurrence-relations/42016/3 "2020-06-25T00:58:05Z")

</div>

Why didn’t I think of that? The actual problem I am solving is this [#227 The Chase - Project Euler](https://projecteuler.net/problem=227)

So wondering if that applies there.

---

<div class="post-metadata">

**Author:** ![dpsanders](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dpsanders/32/3573_2.png) [@dpsanders](https://discourse.julialang.org/u/dpsanders)\
**Post date:** [June 25, 2020, 1:03am UTC](https://discourse.julialang.org/t/what-are-some-good-ways-to-compute-terms-in-recurrence-relations/42016/4 "2020-06-25T01:03:57Z")

</div>

> [@ericphanson](#):
>
> you could just solve for `f(n+1)` and then have a recurrence relation only in terms of lower values,

I don’t think you can do that in this case, since you don’t have enough boundary conditions.

> [@xiaodai](#):
>
> I am just wondering if there’s a good way to solve this in Julia simply. I wrote some recursive function to solve it.

I would be interested to see that function.

What I understand is that you allow positions n=1, \ldots, L, where L is the finish line, and I presume that if you are at n=1 and you try to go left, you stay at n=1 (i.e. you “bounce back” from a wall); in this way there are only a finite number of allowed positions n.

Your f(n) is the _expected_ or _mean_ time to hit the finish line L, starting from site n; this is a _hitting time_ or _(mean) first-passage time_. In particular, f(L) = 0, since you’re already there.

If you write down all the equations for all values of n, you have a system of linear equations for the unknowns f(1), f(2), …, which you can solve using the `\` operator.

A more general technique, in particular, for problems in which there are an infinite number of possible values of n, is to use a _generating function_, i.e. define

F(z) := \sum\_n f(n) z^n

The recurrence relation, together with the boundary conditions, gives you an algebraic equation for F(z), which you can solve. You can then (sometimes) extract the results f(n) as the coefficient of z^n in the resulting expression.

---

<div class="post-metadata">

**Author:** ![xiaodai](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/xiaodai/32/15937_2.png) [@xiaodai](https://discourse.julialang.org/u/xiaodai)\
**Post date:** [June 25, 2020, 1:05am UTC](https://discourse.julialang.org/t/what-are-some-good-ways-to-compute-terms-in-recurrence-relations/42016/5 "2020-06-25T01:05:41Z")

</div>

> [@dpsanders](#):
>
> I would be interested to see that function.

Check it out. It solves the more complicated version not the simplified one I have in my original post.

```julia
# this is an alternate solution that doesn't rely on matrices so is more memory efficient

init_rhs() = begin
    res = Rational{BigInt}[0//1 for i in 1:50]
    res
end

function update227!(rhs_coef, const_part, n, multiplier)
    new_const_part, add_rhs_coef = multiplier .* chase(n)
    rhs_coef .+= add_rhs_coef

    const_part + new_const_part
end

using Memoize

@memoize function chase(dist)
    # C*d_{dist} = A*d_{dist} + B*d_{dist-1} + D*d_{dist-2}
    # lhs_coef = C - A
    lhs_coef = 1//1
    const_part = 1//1

    rhs_coef = init_rhs()

    if dist == 50
        rhs_coef[48] += 2//36
        rhs_coef[49] += 16//36
        rhs_coef[50] += 1//2
    elseif dist == 49
        rhs_coef[47] += 1//36
        rhs_coef[48] += 8//36
        rhs_coef[49] += 1//2 + 1//36
        rhs_coef[50] += 8//36
    elseif dist == 1
        rhs_coef[1] += 1//2 + 1//36
        rhs_coef[2] += 8//36
        rhs_coef[3] += 1//36
    elseif dist == 2
        rhs_coef[1] += 8//36
        rhs_coef[2] += 1//2
        rhs_coef[3] += 8//36
        rhs_coef[4] += 1//36
    else
        rhs_coef[dist-2] += 1//36
        rhs_coef[dist-1] += 8//36
        rhs_coef[dist] += 1//2
        rhs_coef[dist+1] += 8//36
        rhs_coef[dist+2] += 1//36
    end

    if dist != 1
        while any(!=(0), @view rhs_coef[1:dist-1])
            for i in 1:dist-1
                if rhs_coef[i] != 0
                    mult = rhs_coef[i]
                    rhs_coef[i] = 0
                    const_part = update227!(rhs_coef, const_part, i, mult)
                end
            end
        end
    end

    lhs_coef -= rhs_coef[dist]
    rhs_coef[dist] = 0//1

    const_part / lhs_coef, rhs_coef ./ lhs_coef
end

```

---

<div class="post-metadata">

**Author:** ![xiaodai](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/xiaodai/32/15937_2.png) [@xiaodai](https://discourse.julialang.org/u/xiaodai)\
**Post date:** [June 25, 2020, 1:07am UTC](https://discourse.julialang.org/t/what-are-some-good-ways-to-compute-terms-in-recurrence-relations/42016/6 "2020-06-25T01:07:20Z")

</div>

> [@dpsanders](#):
>
> _generating function_ ,

I am well aware of generating functions from my discrete math course. Just wondering if there’s a way to input the recurrence relations coefs and a program can solve it for me using whichever method including generating functions.

---

<div class="post-metadata">

**Author:** ![xiaodai](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/xiaodai/32/15937_2.png) [@xiaodai](https://discourse.julialang.org/u/xiaodai)\
**Post date:** [June 25, 2020, 1:10am UTC](https://discourse.julialang.org/t/what-are-some-good-ways-to-compute-terms-in-recurrence-relations/42016/7 "2020-06-25T01:10:04Z")

</div>

> [@dpsanders](#):
>
> you have a system of linear equations for the unknowns f(1)f(1) , f(2)f(2) , …, which you can solve using the `\` operator.

This was my original idea. However, it fails on large matrices on resource constraint machine (like those on HackerRank) and I just realised that you don’t need a matrix for this more specialised problem, hence the question. You should be able to solve this in a more RAM-efficient way. As I have shown above.

---

<div class="post-metadata">

**Author:** ![dpsanders](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dpsanders/32/3573_2.png) [@dpsanders](https://discourse.julialang.org/u/dpsanders)\
**Post date:** [June 25, 2020, 1:29am UTC](https://discourse.julialang.org/t/what-are-some-good-ways-to-compute-terms-in-recurrence-relations/42016/8 "2020-06-25T01:29:12Z")

</div>

For this particular problem the matrix is symmetric and tridiagonal. Using the `SymTriDiagonal` type should be very efficient (in memory and time) even for large systems.

---

<div class="post-metadata">

**Author:** ![dpsanders](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dpsanders/32/3573_2.png) [@dpsanders](https://discourse.julialang.org/u/dpsanders)\
**Post date:** [June 25, 2020, 1:32am UTC](https://discourse.julialang.org/t/what-are-some-good-ways-to-compute-terms-in-recurrence-relations/42016/9 "2020-06-25T01:32:38Z")

</div>

> [@xiaodai](#):
>
> Just wondering if there’s a way to input the recurrence relations coefs and a program can solve it for me using whichever method including generating functions.

That’s an interesting question (which I don’t have a good answer to). Maybe you could use ModelingToolkit somehow to get the structure of the recurrence relation.

The book [Generatingfunctionology](https://www.math.upenn.edu/~wilf/DownldGF.html) may be useful.

---

<div class="post-metadata">

**Author:** ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Post date:** [June 25, 2020, 8:57am UTC](https://discourse.julialang.org/t/what-are-some-good-ways-to-compute-terms-in-recurrence-relations/42016/10 "2020-06-25T08:57:04Z")

</div>

> [@dpsanders](#):
>
> A more general technique, in particular, for problems in which there are an infinite number of possible values of nn , is to use a _generating function_ , i.e. define

These are very elegant techniques but can run into practical problems when implemented numerically (eg a very nice or even perfectly accurate calculation turning into something less well-conditioned).

Of course these problems can be mitigated, but for small n I would just use the recurrence relation, maybe implemented with a loop, or memoized. Or for non-triangular cases, use a linear system as you suggested.