RogerBW's Blog

I’ve been doing the Weekly Challenges. The latest involved a Fibonacci-like word construction and a search for square-free numbers. (Note that this is open until 6 February 2022.)

Task 1: Fibonacci Words

You are given two strings having same number of digits, $a and $b.

Write a script to generate Fibonacci Words by concatenation of the previous two strings. Finally print 51st digit of the first term having at least 51 digits.

I wondered whether there might be a way of working out the relevant character without actually assembling the strings, but the problem requires the strings to be produced as well, so that's the way I did it. I see no advantage to recursion here. Raku:

sub fibstr($aa,$bb,$limit) {
    my $a=$aa;
    my $b=$bb;
    loop {
        my $c=$a ~ $b;
        say $c;
        if (chars($c) >= $limit) {
            return substr($c,$limit-1,1);
        }
        ($a,$b)=($b,$c);
    }
}

Some of the other languages I'm using distinguish between a 1-length string and a single char value.

Task 2: Square-free Integer

Write a script to generate all square-free integers <= 500.

A square-free integer (or squarefree integer) is defined as an integer which is divisible by no perfect square other than 1. That is, its prime factorization has exactly one factor for each prime that appears in it. For example, 10 = 2 ⋅ 5 is square-free, but 18 = 2 ⋅ 3 ⋅ 3 is not, because 18 is divisible by 9 = 3².

This sequence is in the OEIS, though only up to 113.

The approach I thought of first (there's another one later) is to look at each number and factorise it. In Ruby:

def squarefree(max)

Initialise the output list.

  out=[]

We'll prime-factor each candidate number, with short-cuts. Say I'm considering candidate 500. Since I'm only interested in factors that are square or higher powers of primes, the highest prime that (squared) could be a factor is 19; 23² is 529. So I'll never need a prime larger than the square root of the highest candidate.

  primes=genprimes(Integer.sqrt(max))

And for each individual candidate, I don't need to look at primes higher than its own square root.

  plimit=1

Iterate over each candidate.

  1.upto(max) do |candidate|

A number is square-free if I find no squared prime factors.

    squarefree=true

Increase the limit until the square of the limit is no lower than the candidate.

    while plimit * plimit < candidate do
      plimit += 1
    end

Start factoring here.

    cc=candidate

Check each prime in turn, and break out if we reach the end of the list or exceed the limit.

    primes.each do |pr|
      if pr > plimit then
        break
      end

Count the number of times the number is divisible by this prime.

      n=0
      while cc % pr == 0 do
        n += 1

If there's more than one, it's not square-free, so bail out.

        if n > 1 then
          squarefree=false
          break
        end

Otherwise, reduce the candidate and continue. We don't need to store the prime factors, just confirm that none of them is a square.

        cc=(cc/pr).to_i
      end

If this isn't a square-free number, we don't need to test any more primes.

      if !squarefree then
        break
      end
    end

If we got through all the primes and found nothing to disqualify it, add it to the output list.

    if squarefree then
      out.push(candidate)
    end
  end

Having tested all candidates, return the list of ones that qualified.

  return out
end

Full code on github.

As a bonus for blog readers, I also came up with another algorithm, which starts with an output list containing just 1, then for each prime up to max multiplies each entry in the list by that prime, and stores that in the list if it's not too high. So the list starts as (1); multiplying by 2 makes it (1, 2); multiplying by 3 makes it (1, 2, 3, 6), etc. Because we never multiply by the same prime twice, there can be no non-squarefree numbers in the result.

This doesn't form part of my solution submission, for reasons I'll mention at the end, but it might be useful somewhere. Here's the Rust version, relying on genprimes as above.

fn squarefree(max: u32) -> Vec<u32> {

This time we take primes all the way up to max (because the largest prime within the limit is itself a squarefree number).

    let primes = genprimes(max);

Initialise the output list. The great thing about a BTree is that, while storing to it may be fractionally more expensive than simply slapping a value onto the end of a list, I can cheaply get a list of its members in ascending order. (In fact, I'm already using it in the Rust version of the prime generator.)

    let mut sf = BTreeSet::new();
    sf.insert(1u32);

Loop over each prime, and over each entry in the output list. Note that the list could be changed inside the loop, but in some languages including Rust that will cause this algorithm to go awry; so instead I build a list of values to be inserted, then insert them after the inner loop is complete.

    for pr in primes {
        let mut ii=Vec::new();
        for k in sf.iter() {

Generate a new solution number. If it's low enough, store it (and we'll multiply it in future along with the other entries); if not, bail out, because all the rest of the products in this inner loop will also be too high. (Which is why we need the members of sf in order, because especially with later primes that'll be the majority of the list.)

            let y=k*pr;
            if y <= max {
                ii.push(y);
            } else {
                break;
            }
        }

Insert all the new values into the BTree. (Not that it matters, but none of them will be there already as they're all products of a new prime.)

        if ii.len() > 0 {
            for y in ii {
                sf.insert(y);
            }
        }
    }

Finally, return the BTree as a plain vector of numbers.

    return sf.iter().map(|i| *i).collect::<Vec<u32>>();
}

This gets relatively faster with higher numbers; break-even time seems to be with max somewhere in the range 2^17-2^18, and by about 2^21 on my test machine it's reaching an answer twice as fast as the Rust version of my first solution. However, it does rely for performance on having a working BTree implementation; without that, I'd have to re-sort the keys in sf for each pass through the prime loop, and I don't feel like writing a BTree for the languages that don't have it readily available. (Though I'll probably do it eventually.)