In https://oeis.org/A065069 numbers $n$ such that Fibonacci($n$) is not squarefree, but for all proper divisors $k$ of $n$, Fibonacci($k$) is squarefree, are listed.
OEIS gives a Mathematica program (see link) to calculate these indices. The algorithm used, searches primes $p$ and Fibonacci numbers Fib($k$) such that
$\text{Fib}(k)\;\text{mod}\,p =0$
The numbers $p\cdot k$ are candidates to be indices of primitive squareful Fibonacci numbers.
Can anyone explain me the reasonning behind this algorithm? Probably some theorems about Fibonacci numbers I don't know are used.
I am interested in alternative algorithms as well.