[Blog Post] A Fun Exploration of Perfect, Abundant, and Deficient Numbers

General Usage
, , ,
1

Hey folks,

I recently wrote a blog post with a very small computational treatment of exploring abundant, deficient, and perfect numbers (the latter only briefly) alongside aliquot sequences. This was my first foray into number theory and had a bit of fun with it! Thought I’d share it in case folks find it interesting: A Fun Exploration of Perfect, Abundant, and Deficient Numbers

What was interesting to me was some encounters that I experienced while working with, in particular, aliquot sequences. What I did not realize is that an aliquot sequence can grow to be monstrous in terms of computation. To give an example, I was innocently computing aliquot sequences for several abundant numbers and strangely hit errors when computing on the value n = 138. I didn’t know why I kept getting truncation errors within my computations. I looked up this particular number and discovered that the aliquot sequence for n= 138 has 177 aliquot terms in the sequence! It apparently encounters numbers so big that it was leading to truncation errors.

Furthermore, when tinkering with how to present some of my computational explorations, I decided to use UnicodePlots. This isn’t so much a problem as it is a commendation to say – WOW what a fantastic package! It made showing visualizations terribly easy and fun! Here’s a plot I made just for fun:

Divisor Count for First 1000 Abundant Numbers 
               ┌                                        ┐ 
             1 ┤■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■ 1 000   
             2 ┤■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■ 995     
             3 ┤■■■■■■■■■■■■■■■■■■■■■■ 677                
             4 ┤■■■■■■■■■■■■■■■■■■■■■ 623                 
             5 ┤■■■■■■■■■■ 308                            
             6 ┤■■■■■■■■■■■■■■■■■■■■■■ 672                
             7 ┤■■■■■■■ 216                               
             8 ┤■■■■■■■■■■■ 347                           
             9 ┤■■■■■■■■ 229                              
   Divisors 10 ┤■■■■■■■■■■ 303                            
            11 ┤■■■■ 116                                  
            12 ┤■■■■■■■■■■■ 336                           
            13 ┤■■■ 96                                    
            14 ┤■■■■■■■ 211                               
            15 ┤■■■■■ 139                                 
            16 ┤■■■■■■ 188                                
            17 ┤■■ 68                                     
            18 ┤■■■■■■■ 224                               
            19 ┤■■ 61                                     
            20 ┤■■■■■■■ 201                               
               └                                        ┘ 
                                  Count

Otherwise, like I said in the conclusion of my blog post, I don’t know what other sorts of computational treatments I could do here and how much analogies hold between math domains, but I would love to hear people’s thoughts on that! In terms of Julia specific questions, I was wondering:

  1. How can I handle extremely large numbers within Julia gracefully?
  2. How can I compute upon extremely large numbers within Julia gracefully?

What I mean by huge for example is the 8^{\text{th}} perfect number, 2305843008139952128, and so forth.

Thanks all!

~ tcp :deciduous_tree:

3 Likes
2

Do you mean things like BigInt / BigFloat ?

julia> x = BigInt(2305843008139952128)
2305843008139952128
julia> y = BigFloat(2305843008139952128)
2.305843008139952128e+18
2 Likes
3

Oh yea @sylvaticus – and from there how would I compute on the values? Is it as simple as now proceeding to compute upon them normally?

P.S. I also apologize as I have never actually used these types within Julia before to express values like this.

4

yes, but of course they are much slower than native machine types (e.g. Int64/Float64)

I have recoded the first part of your post using them… I am not a “big” user of them, so I do not think all the specification I added are needed, I think Julia automatically recast operations between Int and BigInt to a BigInt for example, still… this seems to work:

using Pkg
Pkg.activate(temp=true)
Pkg.add("Primes")

import Primes: factor

function divisors(n)
    d = BigInt[BigInt(1)]
    for (p, e) in factor(n)
        t = BigInt[]
        r = BigInt(1)
        for i in BigInt(1):e
            r *= p
            for u in d
                push!(t, u * r)
            end
	    end
	append!(d, t)
    end
    return sort!(d)
end

i = BigInt(1) 
deficient_numbers = BigInt[]
abundant_numbers = BigInt[]
while true
  divisor_sum = divisors(i)[1:end-1] |> sum
  if divisor_sum < i && length(deficient_numbers) != 1000
    push!(deficient_numbers, i)
  elseif divisor_sum > i && length(abundant_numbers) != 1000
    push!(abundant_numbers, i)
  end
  i += BigInt(1)
  length(abundant_numbers) == 1000 && length(deficient_numbers) == 1000 ? break : continue
end
1 Like
5

Oh my goodness! Let me try this with some experiments! Thank you @sylvaticus !

6

You can use the fact that if n = \prod p_i^{e_i} then \sigma(n) = \prod \frac{p_i^{e_i + 1} - 1}{p_i - 1} to avoid allocating a vector to store the divisors. Something like this:

using Primes

function sum_divisors(n)
    s = one(n)
    for (p, e) in Primes.factor(n)
         s *= (p^(e + 1) - 1) ÷ (p - 1)
    end
    s
end

function get_abundant_and_deficient_numbers(n::T) where T <: Integer
    # get the first n abundant and deficient numbers
    n_abundants = 0
    n_deficients = 0

    abundants = sizehint!(T[], n)
    deficients = sizehint!(T[], n)

    k = 1
    while n_abundants < n || n_deficients < n
         σ = sum_divisors(k)
         if σ > 2k && n_abundants < n
             n_abundants += 1
             push!(abundants, k)
        elseif σ < 2k && n_deficients < n
             n_deficients += 1
             push!(deficients, k)
        end

        k += 1
    end

    abundants, deficients
end
1 Like