Problem 193: Squarefree Numbers

Project Euler Problem 193 Solution

EulerSolve provides an optimized solution for Project Euler Problem 193, Squarefree Numbers, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary Let $Q(N)$ denote the number of squarefree positive integers not exceeding $N$. Problem 193 asks for $Q(2^{50})$, where an integer is squarefree if no prime square divides it. The implementations compute this count by turning the condition “has no square prime factor” into an exact divisor sum. This also matches the original wording “below $2^{50}$”, because $2^{50}$ itself is divisible by $4$ and therefore is not squarefree. Mathematical Approach A direct approach would test each integer separately and look for a square divisor. The real solution reverses that viewpoint: instead of scanning numbers one by one, it counts how many integers are multiples of each square and combines those counts with Möbius signs. Detecting square factors with the Möbius function Write $\mu(d)$ for the Möbius function. For a fixed integer $n$, consider $$\sum_{d^2\mid n}\mu(d).$$ Only squarefree values of $d$ can contribute, because $\mu(d)=0$ whenever $d$ contains a repeated prime factor. If $S(n)$ is the set of primes whose squares divide $n$, then every contributing $d$ is the product of a subset of $S(n)$. Therefore $$\sum_{d^2\mid n}\mu(d)=\sum_{T\subseteq S(n)}(-1)^{|T|}=(1-1)^{|S(n)|}.$$ If no prime square divides $n$, then $S(n)$ is empty and the sum is $1$. Otherwise the sum is $0$....

Detailed mathematical approach

Problem Summary

A positive integer $n$ is squarefree if no square of a prime divides $n$. The statement lists $1, 2, 3, 5, 6, 7, 10, 11$ as squarefree and $4, 8, 9, 12$ as not squarefree; $1$ is squarefree because no prime divides it at all. We must find how many squarefree numbers lie below $2^{50}$, that is, among $1,2,\dots,2^{50}-1$. ^[1] Write $Q(N)$ for the number of squarefree $n$ with $1\le n\le N$; the twelve examples of the statement are the numbers $1$ to $12$, so $Q(12)=8$. Since $2^{50}$ is divisible by $4$, it is not squarefree, so the requested count equals $Q(2^{50})$, and from now on $N=2^{50}$.

The whole solution in one sentence: start with all $N$ numbers, remove the multiples of $4$, $9$, $25$ and of every other prime square, put back the numbers that were removed twice, and so on; the signs of these corrections are the values of the Möbius function, so the count becomes a sum with one term for each square $d^2\le N$, only $2^{25}$ terms, and a linear sieve computes all the signs in one pass.

Mathematical Approach

Testing all $2^{50}\approx1.13\cdot10^{15}$ numbers one by one is far too slow. Instead we count, for each square, how many numbers up to $N$ it divides, and combine these counts with the right signs.

1. Small cases

The squarefree numbers up to $30$ are

$1, 2, 3, 5, 6, 7, 10, 11, 13, 14,$
$15, 17, 19, 21, 22, 23, 26, 29, 30,$

so $Q(10)=7$, $Q(12)=8$ and $Q(30)=19$ (OEIS A005117). ^[3] Each of the other eleven numbers, $4, 8, 9, 12, 16, 18,$ $20, 24, 25, 27, 28$, is divisible by $4$, $9$ or $25$. Squares of composite numbers need no separate test: if $d^2$ divides $n$ for some $d\gt1$, then $p^2$ divides $n$ for every prime $p$ that divides $d$. So $n$ is squarefree exactly when no square greater than $1$ divides it, ^[2] and for $n\le N$ only the primes $p\le\sqrt N$ matter, since a larger prime has $p^2\gt N$.

2. Removing the multiples of prime squares

Among $1,2,\dots,N$ the multiples of $m$ are $m,2m,\dots,km$ with $k=\lfloor N/m\rfloor$, where $\lfloor x\rfloor$ is $x$ rounded down. For $N=30$ there are $\lfloor30/4\rfloor=7$ multiples of $4$, $\lfloor30/9\rfloor=3$ multiples of $9$ and $\lfloor30/25\rfloor=1$ multiple of $25$. No number up to $30$ is a multiple of two of these squares, because the smallest such number is $4\cdot9=36$. Removing the three groups therefore removes $11$ different numbers and leaves $30-7-3-1=19$, as in Section 1.

3. Numbers that are removed twice

For $N=100$ the prime squares up to $N$ are $4$, $9$, $25$ and $49$, and removing their multiples leaves $100-25-11-4-2=58$. That is too few: $36$ and $72$ are multiples of both $4$ and $9$, and $100$ is a multiple of both $4$ and $25$, so these three numbers were removed twice. Putting each of them back once gives $58+3=61$, and a direct check confirms $Q(100)=61$.

These numbers can be counted as well. By unique factorisation, $n$ is a multiple of both $p^2$ and $q^2$, for different primes $p$ and $q$, exactly when it is a multiple of $(pq)^2$. Up to $100$ there are $\lfloor100/36\rfloor=2$ multiples of $6^2$ and $\lfloor100/100\rfloor=1$ multiple of $10^2$, while $4\cdot49$ and $9\cdot25$ exceed $100$. In general, $n$ is a multiple of all of $p_1^2,\dots,p_k^2$, for distinct primes $p_i$, exactly when it is a multiple of $(p_1p_2\cdots p_k)^2$.

4. Inclusion–exclusion and the Möbius function

For each prime $p\le\sqrt N$ let $A_p$ be the set of multiples of $p^2$ among $1,\dots,N$. The squarefree numbers are those in none of these sets, and the inclusion–exclusion principle counts them: from all $N$ numbers subtract the sizes of the sets, add the sizes of the intersections of two sets, subtract those of three, and so on. ^[4] By Section 3, the intersection of the sets $A_p$ for the primes of a set $S$ consists of the multiples of $d^2$, where $d$ is the product of these primes; it has $\lfloor N/d^2\rfloor$ elements and enters with the sign $(-1)^{|S|}$. By unique factorisation each product $d$ of distinct primes comes from exactly one set $S$, and the empty set gives $d=1$ and the count $N$.

This sign is the Möbius function: $\mu(d)=(-1)^k$ when $d$ is a product of $k$ distinct primes, so $\mu(1)=1$, and $\mu(d)=0$ when the square of a prime divides $d$. ^[5] Its values for $d=1,\dots,12$ are $1,-1,-1,0,-1,1,-1,0,0,1,-1,0$. Adding the remaining $d$ and dropping those with $d^2\gt N$ changes nothing: $\mu(d)=0$ if $d$ is not a product of distinct primes, and $\lfloor N/d^2\rfloor=0$ if $d^2\gt N$, which holds whenever $d$ has a prime factor above $\sqrt N$. So the sum may run over all $d\le\sqrt N$, and we get the formula used by all three programs: ^[2]

$$Q(N)=\sum_{d=1}^{\lfloor\sqrt N\rfloor}\mu(d)\left\lfloor\frac{N}{d^2}\right\rfloor.$$

5. Every number is counted correctly

The formula can also be checked number by number. A number $n\le N$ is one of the $\lfloor N/d^2\rfloor$ multiples of $d^2$ exactly when $d^2\mid n$, so the sum counts it

$$c(n)=\sum_{d^2\mid n}\mu(d)$$

times. Let $k$ be the number of primes whose square divides $n$. Only the $d$ with $\mu(d)\ne0$ contribute; such a $d$ is a product of distinct primes, and $d^2\mid n$ holds exactly when the squares of all these primes divide $n$. So the contributing $d$ correspond to the subsets of the $k$ primes, a subset of $j$ primes has the sign $(-1)^j$, and by the binomial theorem

$$c(n)=\sum_{j=0}^{k}\binom{k}{j}(-1)^j=(1-1)^k,$$

which is $1$ for $k=0$, that is, for squarefree $n$, and $0$ for $k\ge1$. For example, $72=2^3\cdot3^2$ is counted for $d=1$, removed for $d=2$ and $d=3$ and put back for $d=6$, so $c(72)=1-1-1+1=0$; for $900=2^2\cdot3^2\cdot5^2$ the eight divisors of $30$ give $c(900)=1-3+3-1=0$. Adding $c(n)$ over all $n\le N$ counts each squarefree number once and no other number, which proves the formula a second time.

6. Worked example: $Q(100)$ from the formula

For $N=100$ the sum runs over $d\le10$. With $\mu(1)=1$, $\mu(2)=\mu(3)=\mu(5)=\mu(7)=-1$, $\mu(6)=\mu(10)=1$ and $\mu(4)=\mu(8)=\mu(9)=0$, the terms for $d=1,2,3,5,6,7,10$ give

$$\begin{aligned}Q(100)={}&100-25-11-4\\&+2-2+1=61.\end{aligned}$$

These are the steps of Section 3 in another order: $-2$ removes the multiples of $49$, and $+2$ and $+1$ put back the multiples of $36$ and of $100$. The terms for $d=4$, $8$ and $9$ vanish, because a multiple of $16$, $64$ or $81$ is already a multiple of $4$ or $9$ and has been handled by $d=2$ or $d=3$. In the same way, the $20$ nonzero terms with $d\le32$ give $Q(1024)=624$; the C++ program checks both values before it starts.

7. The density $6/\pi^2$

The formula also shows how large $Q(N)$ is. Write $s=\lfloor\sqrt N\rfloor$ and $\lfloor N/d^2\rfloor=N/d^2-\theta_d$ with $0\le\theta_d\lt1$, where $\theta_1=0$. Then $Q(N)$ differs from $N\sum_{d\le s}\mu(d)/d^2$ by less than $s$. Extending this sum to all $d$ changes it by at most $N\sum_{d\gt s}1/d^2$, which is less than $N/s\le2\sqrt N$, because $1/d^2\lt1/(d-1)-1/d$, these differences add up to $1/s$, and $\sqrt N\lt s+1\le2s$. Multiplying out the product of the factors $1-1/p^2$ over the primes $p\le P$ gives the terms $\mu(d)/d^2$ of all $d$ whose prime factors are at most $P$, as in Section 4; every missing term has $d\gt P$, so the product differs from $\sum_{d\ge1}\mu(d)/d^2$ by at most $\sum_{d\gt P}1/d^2$, which tends to $0$ as $P$ grows. Euler's solution of the Basel problem, $\sum_{n\ge1}1/n^2=\pi^2/6$, together with his product formula for this sum, shows that the product equals $6/\pi^2$. ^[6] Hence

$$\left|Q(N)-\frac{6N}{\pi^2}\right|\lt3\sqrt N,$$

and about $60.79\%$ of all numbers are squarefree. ^[2] For example, $Q(100)=61$ and $Q(10^6)=607926$, while $6N/\pi^2$ is about $60.8$ and $607927.1$.

8. The Möbius values from a linear sieve

For $N=2^{50}$ the formula needs $\mu(d)$ for all $d\le M=2^{25}=33554432$. Instead of factorising each $d$, a sieve fills the whole table in one pass. The linear sieve goes through $i=2,3,\dots,M$ and keeps the list of the primes found so far. ^[7] If no earlier step has marked $i$, then $i$ is prime: it is added to the list and gets $\mu(i)=-1$. Then the sieve takes the primes $p=2,3,5,\dots$ of the list in increasing order, as long as $ip\le M$, and marks $ip$ as composite:

if $p$ does not divide $i$, then $\mu(ip)=-\mu(i)$: the new prime $p$ flips the sign, and if a prime square divides $i$, both values are $0$; this is the rule $\mu(ab)=\mu(a)\mu(b)$ for coprime $a$ and $b$, with $\mu(p)=-1$; ^[5]
if $p$ divides $i$, then $p^2$ divides $ip$, so $\mu(ip)=0$, and the loop over the primes stops for this $i$.

For $M=12$: $i=2$ is prime and marks $4$, with $\mu(4)=0$; $i=3$ is prime and marks $6$, with $\mu(6)=1$, and $9$, with $\mu(9)=0$; $i=4$ marks $8$, with $\mu(8)=0$; $i=5$ is prime, marks $10$, with $\mu(10)=1$, and stops because $15\gt12$; $i=6$ marks $12$, with $\mu(12)=0$; $7$ and $11$ are prime, and nothing else is marked. The values agree with the list in Section 4.

Every composite $c\le M$ is marked exactly once. Let $p$ be its smallest prime factor and $i=c/p$. Then $p\le i$, so $p$ is in the list when $i$ is handled; no prime $q\lt p$ divides $i$, and $iq\lt c\le M$, so the loop for $i$ reaches $p$ and marks $c$. Conversely, if the loop for $i'$ marks $c=i'q$, then no prime below $q$ divides $i'$, or the loop would have stopped there; so $q$ is the smallest prime factor of $c$, and $q=p$, $i'=i$. Since $i\lt c$, the value $\mu(i)$ is known by then; and as every composite $j$ is marked in the loop for $j/r\lt j$, with $r$ its smallest prime factor, a number still unmarked when its turn comes is prime. Each pass of the inner loop marks a new composite or ends the loop, so the inner loop runs fewer than $2M$ times in all, and with the outer loop the sieve takes $O(M)$ time; the sieve of Eratosthenes, which can mark a composite once for each of several prime factors, needs $O(M\log\log M)$. ^[7]

How the Code Works

The C++ program ^[8] keeps its settings in the struct Options: the limit, by default $2^{50}$, and whether to run the checkpoints. parse_arguments accepts --skip-checkpoints and --limit= followed by digits, which parse_u64_after_prefix reads into an unsigned 64-bit number without an overflow check, so a value above $2^{64}-1$ silently wraps around; any other argument, or the limit $0$, ends the program with exit code $1$. Unless they are switched off, run_checkpoints requires solve(100) to be $61$ and solve(1024) to be $624$, the values of Section 6, and otherwise ends the program with exit code $2$. Then main prints the value of solve for the limit.

solve implements Sections 4 and 8. It sets max_k, the bound $M=\lfloor\sqrt N\rfloor$, by converting the limit to long double, taking std::sqrt and truncating to an int, which restricts the option to limits below about $2^{62}$; for $N=2^{50}$ the root is exact, because $2^{50}$ and its root $2^{25}$ are representable and the square root is correctly rounded. ^[9] The sieve uses a vector primes with room reserved for $M/10$ entries, enough for all primes up to $2^{25}$, a byte vector is_composite and a vector mu of signed 8-bit integers; the product $ip$ is formed in 64 bits before it is compared with $M$. The second loop adds $\mu(k)\lfloor N/k^2\rfloor$ for every $k\le M$ with $\mu(k)\ne0$ to the signed 64-bit count. No partial sum exceeds $\sum_k N/k^2=\pi^2N/6\lt2^{51}$ in absolute value, so nothing overflows, and the result is returned as an unsigned 64-bit number.

The Python program has no arguments and no checkpoints: its solve sets limit = 2**50, takes max_k from math.isqrt, which computes the exact integer square root of every nonnegative integer, ^[9] runs the same sieve with the lists mu and primes and the bytearray is_composite, and prints the sum, computed with Python's unbounded integers. The Java program does all the work in main, with the limit fixed at 1L << 50. It takes Math.sqrt of a double, exact here as in C++, and uses an int array for $\mu$, a boolean array comp and an int array of maxK / 5 slots for the primes, filled up to the counter pc. That array is ample for $2^{25}$ but would be too small for small limits such as $100$ or $1024$. The sum is a long, and each product of $\mu$ with the long quotient is formed in 64 bits.

Complexity and Verification

With $M=\lfloor\sqrt N\rfloor$ the sieve takes $O(M)$ steps by Section 8 and the summation one pass over $k\le M$, so the time is $O(\sqrt N)$, and the tables take $O(\sqrt N)$ memory. For $N=2^{50}$ the sieve finds $2063689$ primes, marks each of the $31490742$ composite numbers up to $2^{25}$ once and runs its inner loop $53969081$ times, below the bound $2M$; $20398664$ of the numbers $k\le M$ are squarefree and contribute a term. The C++ vectors take about $80$ MB, two bytes per $k$ and the room reserved for the primes; the Java arrays take about $195$ MB, five bytes per $k$ and the primes; Python's list mu alone holds an 8-byte reference per entry on a 64-bit build. The peak memory measured for the Java and Python processes was about $240$ MB and $400$ MB. In tests on a desktop computer the C++ and Java programs took well under a second and the Python program about six seconds.

Faster methods exist: sieving $\mu$ block by block needs only about $N^{1/4}$ memory, and grouping the large $d$ by the value of $\lfloor N/d^2\rfloor$, with the sums of $\mu$ taken from the Mertens function $\sum_{k\le x}\mu(k)$, needs only about $N^{2/5}$ steps and has counted the squarefree numbers up to $10^{36}$. ^[10]

Besides the two C++ checkpoints, we checked for this article that a direct test of every number agrees with the formula for all $N\le5000$; that the sieve's table agrees with $\mu$ from factorisation for all $d\le200000$; and that the C++ program with --limit= agrees with a direct count for all $N\le2000$ and $200$ random $N\le10^6$, and with the second program below for $15$ larger limits up to $2^{49}$, among them squares and numbers one below a square. All three programs print the same count for $N=2^{50}$, and so do two separately written programs that are not part of the published solutions. One sieves $\mu$ in the manner of Eratosthenes and evaluates the same sum. The other uses no Möbius function: since every $n$ is $a^2b$ with a squarefree $b$ in exactly one way, ^[2] $N=\sum_{a\ge1}Q(\lfloor N/a^2\rfloor)$, and this recursion, with a direct table for arguments up to $2^{30}$, yields $Q(N)$.

Footnotes and References

The references below supply the background facts used above. The inclusion–exclusion formula for $Q(N)$, its check number by number, the bound of Section 7 and the Möbius rule of the linear sieve are derived explicitly in this article.

Project Euler 193 — Squarefree Numbers. The official statement calls a positive integer squarefree if no square of a prime divides it, lists 1, 2, 3, 5, 6, 7, 10, 11 as squarefree and 4, 8, 9, 12 as not, and asks how many squarefree numbers there are below 2⁵⁰.
Wikipedia — Square-free integer. A square-free integer is divisible by no square other than 1, and every positive integer is m²k with a square-free k in exactly one way. The article gives Q(x) = Σ μ(d)⌊x/d²⌋ over d ≤ √x, Q(x) = 6x/π² + O(√x), the density 6/π² ≈ 0.6079 and Q(10⁶) = 607926.
OEIS A005117 — Squarefree numbers. The sequence begins 1, 2, 3, 5, 6, 7, 10, 11, 13, 14, …, and the entry gives the density 6/π².
Wikipedia — Inclusion–exclusion principle. Its complementary form counts the elements of a finite set that lie in none of the subsets A₁, …, Aₙ by an alternating sum of the sizes of their intersections.
Wikipedia — Möbius function. μ(n) is 1 for n = 1, (−1)ᵏ for a product of k distinct primes and 0 when the square of a prime divides n; its first values are 1, −1, −1, 0, −1, 1, −1, 0, 0, 1, −1, 0, and μ(ab) = μ(a)μ(b) for coprime a and b.
Wikipedia — Basel problem. Euler proved that the sum of 1/n² over all positive integers n is π²/6, and the article writes 1 = (π²/6) ∏ₚ (1 − 1/p²) over all primes p, the Euler product for 1/ζ(2).
cp-algorithms — Linear Sieve. Each number i is multiplied by the primes up to its smallest prime factor, so every composite number is reached exactly once and the sieve runs in O(n) time, against O(n log log n) for the sieve of Eratosthenes; the page attributes the algorithm to Paul Pritchard.
C++ — Euler193.cpp. The linked immutable C++ revision contains parse_arguments, solve and the checkpoints Q(100) = 61 and Q(1024) = 624 of run_checkpoints. It is the implementation record for this article.
cppreference — std::sqrt: correctly rounded, as IEEE 754 requires, so the result is exact whenever it is representable; Java — Math.sqrt: the double value closest to the true square root; Python — math.isqrt: the integer square root of a nonnegative integer.
Jakub Pawlewicz — Counting Square-Free Numbers (arXiv:1107.4890). Starting from the same formula, the paper evaluates the direct sum in Õ(√n) time, where Õ ignores logarithmic factors, with O(n^(1/4)) memory by sieving μ in blocks, and in Õ(n^(2/5)) time by grouping the large d by the value of ⌊n/d²⌋ and using the Mertens function; it reports the count up to 10³⁶.

Problem 193 source code

C++

#include <cmath>
#include <cstdint>
#include <iostream>
#include <string>
#include <vector>

namespace {

using u64 = std::uint64_t;
using i64 = std::int64_t;

struct Options {
    u64 limit = (1ULL << 50);
    bool run_checkpoints = true;
};

bool parse_u64_after_prefix(const std::string& arg, const std::string& prefix, u64& value) {
    if (arg.rfind(prefix, 0U) != 0U) {
        return false;
    }
    const std::string tail = arg.substr(prefix.size());
    if (tail.empty()) {
        return false;
    }

    u64 parsed = 0;
    for (char c : tail) {
        if (c < '0' || c > '9') {
            return false;
        }
        parsed = parsed * 10ULL + static_cast<u64>(c - '0');
    }
    value = parsed;
    return true;
}

bool parse_arguments(int argc, char** argv, Options& options) {
    for (int i = 1; i < argc; ++i) {
        const std::string arg(argv[i]);

        if (arg == "--skip-checkpoints") {
            options.run_checkpoints = false;
            continue;
        }
        if (parse_u64_after_prefix(arg, "--limit=", options.limit)) {
            continue;
        }

        std::cerr << "Unknown argument: " << arg << '\n';
        return false;
    }

    return options.limit >= 1;
}

u64 solve(const u64 limit) {
    const int max_k = static_cast<int>(std::sqrt(static_cast<long double>(limit)));

    std::vector<int> primes;
    primes.reserve(static_cast<std::size_t>(max_k / 10));
    std::vector<std::uint8_t> is_composite(static_cast<std::size_t>(max_k + 1), 0);
    std::vector<std::int8_t> mu(static_cast<std::size_t>(max_k + 1), 0);

    mu[1] = 1;
    for (int i = 2; i <= max_k; ++i) {
        if (!is_composite[static_cast<std::size_t>(i)]) {
            primes.push_back(i);
            mu[static_cast<std::size_t>(i)] = -1;
        }

        for (int p : primes) {
            const i64 v = static_cast<i64>(i) * p;
            if (v > max_k) {
                break;
            }
            is_composite[static_cast<std::size_t>(v)] = 1;
            if ((i % p) == 0) {
                mu[static_cast<std::size_t>(v)] = 0;
                break;
            }
            mu[static_cast<std::size_t>(v)] = static_cast<std::int8_t>(-mu[static_cast<std::size_t>(i)]);
        }
    }

    i64 count = 0;
    for (i64 k = 1; k <= max_k; ++k) {
        const std::int8_t m = mu[static_cast<std::size_t>(k)];
        if (m == 0) {
            continue;
        }
        count += static_cast<i64>(m) * static_cast<i64>(limit / static_cast<u64>(k * k));
    }

    return static_cast<u64>(count);
}

bool run_checkpoints() {
    if (solve(100ULL) != 61ULL) {
        std::cerr << "Checkpoint failed for limit=100" << '\n';
        return false;
    }
    if (solve(1024ULL) != 624ULL) {
        std::cerr << "Checkpoint failed for limit=1024" << '\n';
        return false;
    }
    return true;
}

}  // namespace

int main(int argc, char** argv) {
    Options options;
    if (!parse_arguments(argc, argv, options)) {
        return 1;
    }
    if (options.run_checkpoints && !run_checkpoints()) {
        return 2;
    }

    std::cout << solve(options.limit) << '\n';
    return 0;
}

Python

# Problem 193: Squarefree Numbers
# Count squarefree numbers below 2^50 using Mobius function.

from math import isqrt

def solve():
    limit = 2**50
    max_k = isqrt(limit)
    # Linear sieve for Mobius function
    mu = [0]*(max_k+1)
    mu[1] = 1
    primes = []
    is_composite = bytearray(max_k+1)
    for i in range(2, max_k+1):
        if not is_composite[i]:
            primes.append(i)
            mu[i] = -1
        for p in primes:
            if i*p > max_k: break
            is_composite[i*p] = 1
            if i % p == 0:
                mu[i*p] = 0
                break
            mu[i*p] = -mu[i]
    count = 0
    for k in range(1, max_k+1):
        if mu[k] != 0:
            count += mu[k] * (limit // (k*k))
    print(count)

solve()

Java

public class Euler193 {
    public static void main(String[] args) {
        long limit = 1L << 50;
        int maxK = (int) Math.sqrt((double) limit);
        int[] mu = new int[maxK + 1];
        mu[1] = 1;
        boolean[] comp = new boolean[maxK + 1];
        int[] primes = new int[maxK / 5];
        int pc = 0;
        for (int i = 2; i <= maxK; i++) {
            if (!comp[i]) {
                primes[pc++] = i;
                mu[i] = -1;
            }
            for (int j = 0; j < pc; j++) {
                long v = (long) i * primes[j];
                if (v > maxK)
                    break;
                comp[(int) v] = true;
                if (i % primes[j] == 0) {
                    mu[(int) v] = 0;
                    break;
                }
                mu[(int) v] = -mu[i];
            }
        }
        long count = 0;
        for (long k = 1; k <= maxK; k++)
            if (mu[(int) k] != 0)
                count += mu[(int) k] * (limit / (k * k));
        System.out.println(count);
    }
}