I'm not an expert of number-theoretic algorithms, but it seems to me that you can employing the Chinese Remainder Theorem to obtain a decent first stab at a "sieving" algorithm.
Use several registers r[p], each storing a residue modulo p2 for some prime p. Define such a register for various small primes p ∈ {2, 3, 5, ... , P} for some suitably large prime P. You will use these to represent an integer R, such that R ≡ r[p] (mod p2). You shouldn't need too many such registers to faithfully represent even reasonably large non-negative numbers (the registers can uniquely determine any integer from 1 to 22 × 32 × 52 × ... × P2).
Whenever you wish to increment R, increment each of the registers r[p] as well. If R is square-free, none of these registers will be zero modulo the square of the appropriate prime. For sufficiently small integers N, you can even characterize N as being square-free if none of these residues are zero. Put another way, the more registers r[p] you maintain, the larger the range of square-free numbers you can automatically detect using these registers.
Suppose that you want to test for square-freeness up to some upper limit N. What do you need to test square-freeness using nothing but registers such as I've described? What you need is for any composite number less than N to have a prime factor from the list of registers that you maintain; that is, you need registers for each of the primes up to √N.
If you're going to iterate through all integers anyway, you can discover the list of primes that you need to characterize square-freeness at the same time by the Sieve of Erasthotenes, and construct the list of residues modulo squares-of-primes as you go. Any time you find a new prime P, so long as P2 < N, add P to the list of primes for which you maintain registers and initialize a register for residues modulo P2.
As Ross suggests in his answer, you can use results on the distribution of square-free numbers to obtain an upper bound for the Nth square-free number (it should grow more slowly than N ln(N) or so, in any case). This gives you an upper bound on the number of registers that you have to maintain as you search for square-free integers.
From this, you should actually be able to show a polynomial-time asymptotic upper bound on the runtime of your procedure.