Solution 1:

The number of divisors of a square is the divisor function convolved with the square of the Möbius function (see $g(n)$ here)

$$e(n)=d(n^2)=(d\star\mu^2)(n)$$

by combining these three identities from the linked table

$$e=1\star1\star Q_2,\space d=1\star1,\space Q_2=\mu^2$$

and since

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

therefore

$$e(n)=\sum_{a \left| n \right.} d \left( \frac{n}{a} \right) \sum_{b^2 \left| a \right.} \mu \left( b \right)$$

which can be simplified and rewritten as

$$e(n)=\sum_{b^2 \left| n \right.} \mu \left( b \right) d_3 \left( \frac{n}{b^2} \right)$$

where $d_3(n)$ is the number of ways that a given number can be written as a product of three integers. This identity can be verified by noting that $e(n)$ is multiplicative and checking at prime powers which yields $e(p^a)={2a+1}$ and can be compared with $d(p^a)={a+1}$. In particular note that $d_3(p^a)={\binom{a+2}{2}}$ (see $d_k$ here).

Then the summation of the number of divisors of the square numbers can be computed as:

$$E(x)=\sum_{n\le x} d(n^2) =\sum_{n \leq x} \sum_{b^2 \left| n \right.} \mu \left( b \right) d_3 \left( \frac{n}{b^2} \right)$$

which can be reorganized as:

$$E(x)=\sum_{b \leq \sqrt{x}} \mu \left( b \right) \sum_{n \leq x / b^2} d_3 \left( n \right)$$ $$E(x)=\sum_{a \leq \sqrt{x}} \mu \left( a \right) D_3 \left( \frac{x}{a^2} \right)$$

where $D_3$ is the summatory function for $d_3$. Since $D_3(x)$ can be computed in $O(x^{2/3})$ time complexity using the three dimensional analogue of the hyperbola method, $E(x)$ can also therefore be computed in

$$O(\int_{a=1}^\sqrt{x}{(\frac{x}{a^2})}^{2/3}da)=O(x^{2/3})$$

which is better than $O(x)$ as desired.

By taking an $O(x^{1/3})$ algorithm to compute $D(x)$ and using it for $D_3(x)$, this bound can be reduced to $O(x^{5/9})$. You can find such an algorithm and one formula for $D_3(x)$ in my article here which uses the notation $T(n)$ and $T_3(n)$.