Goodstein Sequence until 50th term

goodstein_sequence.jl

# Ref.: https://en.wikipedia.org/wiki/Goodstein%27s_theorem#Goodstein_sequences

@kwdef struct GoodsteinTerm
    coefficient::Int
    base::Int
    exp::Int
    value::Function = () -> coefficient * base ^ exp
end

function factorize(n::Int, base::Int, tmp_arr=[])
  if n < base
    return [tmp_arr..., n % base]
  end
  return factorize(Int(floor(n / base)), base, [tmp_arr..., n % base])
end

function factors_decomposition(n::Int, base::Int)::Vector{GoodsteinTerm}
    factors = factorize(n, base)

    gen_terms(acc, (idx, val)) = val != 0 ? append!(acc, [GoodsteinTerm(coefficient=val, base=base, exp=(idx - 1))]) : acc

    reduce(gen_terms, enumerate(factors), init = [])
end

b₀ = 2
n = 4
tmp_factors = factors_decomposition(n, b₀)

# Just calculate the first 50th term, because this sequence grows too fast and take too long to stop at 0 
while n > 0 && b₀ <= 51
    # Check if the exponent's value is equal the base
    new_exp(curr, base) = curr == base ? curr + 1 : curr

    # Increment the base (if necessary the exponent) and conserve the coefficients
    new_factor(curr) = GoodsteinTerm(coefficient = curr.coefficient,
                                     base = curr.base + 1,
                                     exp = new_exp(curr.exp, b₀))

    global tmp_factors = [new_factor(factor) for factor in tmp_factors]

    global b₀ += 1

    new_value = reduce((acc, curr) -> curr.value() + acc, tmp_factors, init=0) - 1

    global tmp_factors = factors_decomposition(new_value, b₀)
    global n = new_value

    println("$(b₀ - 2) => $(new_value)")
end