# Trying Julia for Analytic combinatorics

**URL:** <https://discourse.julialang.org/t/trying-julia-for-analytic-combinatorics/103682>\
**Category:** General Usage\
**Tags:** symbolics, probability, dynamical-systems, inf\
**Created:** [September 8, 2023, 9:52pm UTC](https://discourse.julialang.org/t/trying-julia-for-analytic-combinatorics/103682 "2023-09-08T21:52:07Z")\
**Posts on this page:** 11\
**Page:** 1

<div class="post-metadata">

**Author:** ![fargolo](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fargolo/32/37515_2.png) [@fargolo](https://discourse.julialang.org/u/fargolo)\
**Post date:** [September 8, 2023, 9:52pm UTC](https://discourse.julialang.org/t/trying-julia-for-analytic-combinatorics/103682/1 "2023-09-08T21:52:07Z")

</div>

Hello, fellow colleagues in the Julia community.

This question may not be a good fit for the Julia Discourse, but any help is appreciated.

I was exploring Analytic Combinatorics (by Flajolet & Sedgewick). The part about binary words (words over a binary alphabet). contains the following generating function. in eq. 54 (pg. 52)

![image](https://global.discourse-cdn.com/julialang/original/3X/a/9/a975504e0e60b9d53cee655d1339a576dc0277ac.png)

The following part in page 53 (equation without number) describes the probabilities associated with the maximum length of runs of consecutive letters:

 ![image](https://global.discourse-cdn.com/julialang/original/3X/0/f/0f6edc52845880db5dd72ec987c0f4002b74b065.png)

I used SymPy.jl to implement it:

```julia
using SymPy
#using FastRationals

@vars z r n

# Binary words that never have more than r consecutive identical letters is found to be (set α = β = r)
# Flajolet pag. 52
w_rr(r,z) = (1-z^(r+1))/(1-2z+z^(r+1)) # OGF
w_rr(r,z) = sum(z^x for x in 0:r)/(1 - sum(z^x for x in 1:r)) # Alternate form

function n_words(r)
    coefs = collect(series(w_rr(r,z),z,0,r+1),z)
    coefs.coeff(z,r)
end

function p_k(k,n)
    #a = FastRational{Int128}(1/(2^n))
    #a = Rational{BigInt}(1/(2^n))
    a = 1/(2^n)
    result_sym = a*(n_words(k) - n_words(k-1))
    return(result_sym)
end

```

Testing for `k=6` and verifying the Taylor series for `r=15`:

```julia
julia> p_k(6,n)
    -n
32⋅2
julia>series(w_rr(15,z),z,0,15+1)
             2 3 4 5 6 7 8 9 10 11 12 13 14        
1 + 2⋅z + 4⋅z + 8⋅z + 16⋅z + 32⋅z + 64⋅z + 128⋅z + 256⋅z + 512⋅z + 1024⋅z + 2048⋅z + 4096⋅z + 8192⋅z + 16384⋅z + 32768

  15 ⎛ 16⎞
⋅z + O⎝z ⎠

```

My problems:

1 - Running `p_k(6,200)` directly returns ∞ or another undefined value, even in with `FastRational{Int128}(1/(2^n))` or `Rational{BigInt}(1/(2^n))` that are commented out. I can run `p_k(6,n)` passing symbolic `n` and then evaluate the output directly in the REPL (`julia>32*2^(-200)`) though.

2 - The result seems to be in accordance with the formula `(1/2^n)*([z^6]W<6,6> - [z^5]W<5,5>)`, but it is quite different from ‘template’ for k=6. The coefs. in the Taylor expansion seem to be just `[z^n]W<k,k> = 2^n`. Something must be wrong, since monotonically increasing values for k would not give the `P(K)` distribution presented in the table for integers k in [1,12].

3 - I intend to calculate probabilities and test hypotheses about structures in a [recurrence matrix](https://juliadynamics.github.io/DynamicalSystems.jl/v1.1/rqa/rplots/) (trapping time, vertical, diagonal), which are also binary words. Any thoughts on that?

---

<div class="post-metadata">

**Author:** ![Dan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dan/32/42581_2.png) [@Dan](https://discourse.julialang.org/u/Dan)\
**Post date:** [September 8, 2023, 11:09pm UTC](https://discourse.julialang.org/t/trying-julia-for-analytic-combinatorics/103682/2 "2023-09-08T23:09:59Z")

</div>

Maybe `Rational{BigInt}(1/(2^n))` doesn’t give the expected result, as the `2` in `2^n` is parsed as an Int64 and when raised to 200 gives 0 (i.e. `zero(Int)`) making the whole expression `1//0`. The intended value can be obtained with `1/BigInt(2)^200`:

```julia
julia> Rational{BigInt}(1/(2^200)) # bad result
1//0

julia> 1//(BigInt(2)^200) # good result
1//1606938044258990275541962092341162602522202993782792835301376

```

---

<div class="post-metadata">

**Author:** ![fargolo](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fargolo/32/37515_2.png) [@fargolo](https://discourse.julialang.org/u/fargolo)\
**Post date:** [September 8, 2023, 11:27pm UTC](https://discourse.julialang.org/t/trying-julia-for-analytic-combinatorics/103682/3 "2023-09-08T23:27:17Z")

</div>

Thank you!  
It worked like a charm.

Do you have any idea about my 2nd and 3rd questions?

---

<div class="post-metadata">

**Author:** ![fargolo](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fargolo/32/37515_2.png) [@fargolo](https://discourse.julialang.org/u/fargolo)\
**Post date:** [September 9, 2023, 2:30am UTC](https://discourse.julialang.org/t/trying-julia-for-analytic-combinatorics/103682/4 "2023-09-09T02:30:55Z")

</div>

I fixed it thanks to [RicBit](https://github.com/ricbit).

`W_n<k,k>` notation stands for `[z^n] Taylor(W<k,k>)`. That is, we need the 200th coefficient of the Taylor expansion of W\<k,k\>.

This version works:

```julia
function n_words(r;n_tot=200)
    coefs = collect(series(w_rr(r,z),z,0,n_tot+1),z)
    coefs.coeff(z,n_tot)
end

function p_k(k,n)
    #a = FastRational{Int128}(1/(2^n))
    #a = Rational{BigInt}(1/(2^n))
    a = 1//(BigInt(2)^n)
    #a = 1/(2^n)
    result_sym = a*(n_words(k) - n_words(k-1))
    return(Float64(result_sym))
end

julia>for i in 3:12
    print(p_k(i,200))
    print("\n")
end

# Correct results
6.549422864002916e-8
0.0007075224314585063
0.033979610958548824
0.1660282049825555
0.2574174128392719
0.22352195216110582
0.14594875198080856
0.08292457405578264
0.0440234813776148
0.022607888203734293

```

---

<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:** [September 9, 2023, 6:35am UTC](https://discourse.julialang.org/t/trying-julia-for-analytic-combinatorics/103682/5 "2023-09-09T06:35:35Z")

</div>

This is nice, but I don’t think using SymPy is really using Julia 😉

Maybe [GitHub - JuliaDiff/TaylorSeries.jl: Taylor polynomial expansions in one and several independent variables.](https://github.com/JuliaDiff/TaylorSeries.jl) can help. (Disclaimer: I am one of the co-authors.)

---

<div class="post-metadata">

**Author:** ![fargolo](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fargolo/32/37515_2.png) [@fargolo](https://discourse.julialang.org/u/fargolo)\
**Post date:** [September 9, 2023, 10:39pm UTC](https://discourse.julialang.org/t/trying-julia-for-analytic-combinatorics/103682/6 "2023-09-09T22:39:51Z")

</div>

I gave it a shot, but it seems that there’s something wrong with the types again.  
The results are blowing up.

```julia
using TaylorSeries

w_rr(r,z) = (1-z^(r+1))/(1-2z+z^(r+1))
n_words(r;n_tot=200) = getcoeff(taylor_expand(z -> w_rr(z,r),order=n_tot+1),n_tot)

function p_k(k,n)
    #a = FastRational{Int128}(1/(2^n))
    #a = Rational{BigInt}(1/(2^n))
    a = 1//(BigInt(2)^n)
    #a = 1/(2^n)
    result_sym = a*(n_words(k;n_tot=n) - n_words(k-1;n_tot=n))
    return(Float64(result_sym))
end

julia>for i in 3:12
    print(p_k(i,200))
    print("\n")
end

3.1854727641374944e6
5.919210583485605e18
2.9142600879986195e27
1.2941726894119323e34
2.808302147070106e39
7.566803340725184e43
4.4580962403509384e47
8.293787625532763e50
6.26377090076716e53
2.29651524957329e56

```

---

<div class="post-metadata">

**Author:** ![fargolo](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fargolo/32/37515_2.png) [@fargolo](https://discourse.julialang.org/u/fargolo)\
**Post date:** [March 12, 2026, 1:15pm UTC](https://discourse.julialang.org/t/trying-julia-for-analytic-combinatorics/103682/7 "2026-03-12T13:15:33Z")

</div>

Late update:  
This turned into a full package for analytic combinatorics

> **[GitHub - fargolo/AnalyticComb.jl: Solutions for combinatorial problems using...](https://github.com/fargolo/AnalyticComb.jl)**
>
> Solutions for combinatorial problems using symbolic methods.

and applications for time-series.

> **[GitHub - fargolo/SymbolicInference.jl: Probabilistic inference using the symbolic method.](https://github.com/fargolo/SymbolicInference.jl)**
>
> Probabilistic inference using the symbolic method.

@dpsanders I eventually incorporated TaylorSeries.jl and it gave a massive speed-up that enabled SymbolicInference to handle the data-intensive time-series analysis using recurrence matrices.

---

<div class="post-metadata">

**Author:** ![Datseris](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/datseris/32/13406_2.png) [@Datseris](https://discourse.julialang.org/u/Datseris)\
**Post date:** [March 12, 2026, 1:28pm UTC](https://discourse.julialang.org/t/trying-julia-for-analytic-combinatorics/103682/8 "2026-03-12T13:28:58Z")

</div>

Interesting! WOuld you be interested in giving a talk at an upcoming [JuliaDynamics monthly meeting](https://discourse.julialang.org/t/juliadynamics-monthly-meetings-round-2/132357) about this?

---

<div class="post-metadata">

**Author:** ![fargolo](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fargolo/32/37515_2.png) [@fargolo](https://discourse.julialang.org/u/fargolo)\
**Post date:** [March 12, 2026, 4:26pm UTC](https://discourse.julialang.org/t/trying-julia-for-analytic-combinatorics/103682/9 "2026-03-12T16:26:39Z")

</div>

Yes, @Datseris !

I will take a look at the link. Thank you.

---

<div class="post-metadata">

**Author:** ![fargolo](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fargolo/32/37515_2.png) [@fargolo](https://discourse.julialang.org/u/fargolo)\
**Post date:** [May 11, 2026, 8:58pm UTC](https://discourse.julialang.org/t/trying-julia-for-analytic-combinatorics/103682/10 "2026-05-11T20:58:51Z")

</div>

Another update:

> **[Combinatorial time-loops: probabilistic inference on time-series based on...](https://link.springer.com/article/10.1140/epjs/s11734-026-02334-7)**
>
> Recurrence quantification analysis (RQA) is inspired by Poincaré’s early studies and describes non-linear characteristics of dynamical systems. It achieves this by identifying similarities between states, pairing each observation with every other....

A white-paper outlining this framework was published in an issue of the \* [The European Physical Journal](https://link.springer.com/journal/11734)

---

<div class="post-metadata">

**Author:** ![fargolo](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fargolo/32/37515_2.png) [@fargolo](https://discourse.julialang.org/u/fargolo)\
**Post date:** [July 9, 2026, 10:22pm UTC](https://discourse.julialang.org/t/trying-julia-for-analytic-combinatorics/103682/11 "2026-07-09T22:22:11Z")

</div>

The JuliaDynamics talk is tomorrow 2 PM London time at:

> **[Join conversation](https://teams.microsoft.com/dl/launcher/launcher.html?url=%2F_%23%2Fl%2Fmeetup-join%2F19%3Ameeting_NGNkMTAxMzAtYTFmYi00NTFiLTg3YmMtNTBiNWE2NzI3MzY1%40thread.v2%2F0%3Fcontext%3D%257b%2522Tid%2522%253a%2522912a5d77-fb98-4eee-af32-1334d8f04a53%2522%252c%2522Oid%2522%253a%25220e929467-a486-4dd6-aadf-a47fca854b85%2522%257d%26anon%3Dtrue&type=meetup-join&deeplinkId=9d04971c-d1d4-4d29-804d-c3093c4a45aa&directDl=true&msLaunch=true&enableMobilePage=true)**