2 deleted 10 characters in body

This is really a response to Dror Speiser's comment (but it is too long to give as a comment), and gives an analytic $O(n^{1/3+\epsilon})$ time and $O(n^{1/6+\epsilon})$ space algorithm to count the number of square free numbers less than $x$ (thus improving on the complexity of the algorithm given in the Jakub Pawlewicz paper cited in the question). If we do exactly in the same way as in my other answer for $$\sum_{n\leq x} |\mu(n)|$$ instead of $$\sum_{n\leq x} \mu(n)$$ by using the fact that $$\frac{\zeta(s)}{\zeta(2s)} = \sum_{n=1}^\infty |\mu(n)| n^{-s}$$ it is true that the time complexity will be $O(x^{1/2+\epsilon})$. There is however a way to combine the analytic method (which is essentially the same as in Lagarias-Odlyzkos 1987 paper "Computing $\pi(x)$ an analytic method" http://www.dtc.umn.edu/~odlyzko/doc/arch/analytic.pi.of.x.pdf) with the elementary identity $$S(n)=\sum_{d=1}^{\lfloor \sqrt n\rfloor} \mu(d) \lfloor \frac n {d^2} \rfloor$$ that is Pawleviec's starting point in obtaining his algorithm. Let $\Phi$ be defined as in my other answer, i.e. let $\Phi(x)$ be a smooth test function such that $\Phi(x)=1$ for $x<0$ and $\Phi(x)=0$ for $x>1$. Consider the identity $$\sum_{n\leq x} |\mu(n)| =\sum_{n=1}^\infty |\mu(n)| \Phi \left( \frac {n-x} {x^{2/3}} \right ) - \sum_{x< n < x + x^{2/3} } |\mu(n)| \Phi \left(\frac {n-x} {x^{2/3}} \right)$$

The first sum can be estimated by a complex integral of length $x^{1/3+\epsilon}$ and by the Odlyzko-SchÃ¶nhage algorithm it can be calculated in $O(x^{1/3+\epsilon})$ time (and $O(x^{1/6+\epsilon})$ space). The remaining sum will have length $O(x^{2/3})$ so might at first glance not be so simple to treat. I claim however that in fact it can be calculated in $O(x^{1/3+\epsilon})$ time. This is where using the elementary identity is handy. Let $\Psi(x)=0$ for $x<0$ and $\Psi(x)=\Phi(x)$ for $x \geq 0$ (this will have a discontinuity at $x=0$) By a similar elementary identity we obtain

$$\sum_{x< n < x + x^{2/3} } |\mu(n)| \Phi \left(\frac {n-x} {x^{2/3}} \right) = \sum_{d=1}^{\lfloor \sqrt {x+x^{2/3}}\rfloor} \mu(d) \sum_{k=1}^\infty \Psi \left(\frac {d^2k-x} {x^{2/3}} \right).$$ For any $d$ the inner sum can be calculated fast (certainly in $O(x^\epsilon)$ time), by either noticing that the sum is empty, just contains one element or the Poisson summation formula. The trick now is to see that there will only be $O(x^{1/3})$ integers $d$ where the inner sum (in $k$) is non empty, namely $d$ must either be $O(x^{1/3})$ or belong to an interval $(\sqrt{ x/k},\sqrt{(x+x^{2/3})/k})$ for $k \leq x^{1/3}$ and there are only $O(x^{1/3})$ such $d$. Now if we can calculate the values of $\mu(d)$ fast (simple sieving seems more difficult in this case) we obtain the desired algorithm. Determining $\mu(d)$ basically comes down to factoring the number. However it is sufficient to determine the factors less than $d^{1/3}$ since if we know that all prime factors are greater than $d^{1/3}$, then either it has one or two prime factors. We can determine if it is a prime fast by the Agrawal-Kayal-Saxena Algorithm algorithm, then of course $\mu(d)=-1$. If it has two prime factors either $\mu(d)=0$, but then it has to be a square, and it can be easily checked, or $\mu(d)=1$. Now factoring the $O(x^{1/3})$ numbers of size $O(x^{1/2})$ into all prime factors less than $x^{1/6}$ will take $O(x^{1/3+\epsilon})$ time and $O(x^{1/6+\epsilon})$ space by Bernstein's algorithm http://cr.yp.to/papers/sf.pdf. I might possibly write up this answer as a paper proper. If anyone has seen these arguments anywhere else or have references I would be interested.

1

This is really a response to Dror Speiser's comment (but it is too long to give as a comment), and gives an analytic $O(n^{1/3+\epsilon})$ time and $O(n^{1/6+\epsilon})$ space algorithm to count the number of square free numbers less than $x$ (thus improving on the complexity of the algorithm given in the Jakub Pawlewicz paper cited in the question). If we do exactly in the same way as in my other answer for $$\sum_{n\leq x} |\mu(n)|$$ instead of $$\sum_{n\leq x} \mu(n)$$ by using the fact that $$\frac{\zeta(s)}{\zeta(2s)} = \sum_{n=1}^\infty |\mu(n)| n^{-s}$$ it is true that the time complexity will be $O(x^{1/2+\epsilon})$. There is however a way to combine the analytic method (which is essentially the same as in Lagarias-Odlyzkos 1987 paper "Computing $\pi(x)$ an analytic method" http://www.dtc.umn.edu/~odlyzko/doc/arch/analytic.pi.of.x.pdf) with the elementary identity $$S(n)=\sum_{d=1}^{\lfloor \sqrt n\rfloor} \mu(d) \lfloor \frac n {d^2} \rfloor$$ that is Pawleviec's starting point in obtaining his algorithm. Let $\Phi$ be defined as in my other answer, i.e. let $\Phi(x)$ be a smooth test function such that $\Phi(x)=1$ for $x<0$ and $\Phi(x)=0$ for $x>1$. Consider the identity $$\sum_{n\leq x} |\mu(n)| =\sum_{n=1}^\infty |\mu(n)| \Phi \left( \frac {n-x} {x^{2/3}} \right ) - \sum_{x< n < x + x^{2/3} } |\mu(n)| \Phi \left(\frac {n-x} {x^{2/3}} \right)$$

The first sum can be estimated by a complex integral of length $x^{1/3+\epsilon}$ and by the Odlyzko-SchÃ¶nhage algorithm it can be calculated in $O(x^{1/3+\epsilon})$ time (and $O(x^{1/6+\epsilon})$ space). The remaining sum will have length $O(x^{2/3})$ so might at first glance not be so simple to treat. I claim however that in fact it can be calculated in $O(x^{1/3+\epsilon})$ time. This is where using the elementary identity is handy. Let $\Psi(x)=0$ for $x<0$ and $\Psi(x)=\Phi(x)$ for $x \geq 0$ (this will have a discontinuity at $x=0$) By a similar elementary identity we obtain

$$\sum_{x< n < x + x^{2/3} } |\mu(n)| \Phi \left(\frac {n-x} {x^{2/3}} \right) = \sum_{d=1}^{\lfloor \sqrt {x+x^{2/3}}\rfloor} \mu(d) \sum_{k=1}^\infty \Psi \left(\frac {d^2k-x} {x^{2/3}} \right).$$ For any $d$ the inner sum can be calculated fast (certainly in $O(x^\epsilon)$ time), by either noticing that the sum is empty, just contains one element or the Poisson summation formula. The trick now is to see that there will only be $O(x^{1/3})$ integers $d$ where the inner sum (in $k$) is non empty, namely $d$ must either be $O(x^{1/3})$ or belong to an interval $(\sqrt{ x/k},\sqrt{(x+x^{2/3})/k})$ for $k \leq x^{1/3}$ and there are only $O(x^{1/3})$ such $d$. Now if we can calculate the values of $\mu(d)$ fast (simple sieving seems more difficult in this case) we obtain the desired algorithm. Determining $\mu(d)$ basically comes down to factoring the number. However it is sufficient to determine the factors less than $d^{1/3}$ since if we know that all prime factors are greater than $d^{1/3}$, then either it has one or two prime factors. We can determine if it is a prime fast by the Agrawal-Kayal-Saxena Algorithm algorithm, then of course $\mu(d)=-1$. If it has two prime factors either $\mu(d)=0$, but then it has to be a square, and it can be easily checked, or $\mu(d)=1$. Now factoring the $O(x^{1/3})$ numbers of size $O(x^{1/2})$ into all prime factors less than $x^{1/6}$ will take $O(x^{1/3+\epsilon})$ time and $O(x^{1/6+\epsilon})$ space by Bernstein's algorithm http://cr.yp.to/papers/sf.pdf. I might possibly write up this answer as a paper proper. If anyone has seen these arguments anywhere else or have references I would be interested.