Determine all functions $f: \mathbb{R} \rightarrow \mathbb{R}$ that satisfy 
$$ (f(a)-f(b))(f(b)-f(c))(f(c)-f(a))=f\left(a b^{2}+b c^{2}+c a^{2}\right)-f\left(a^{2} b+b^{2} c+c^{2} a\right) $$ 
for all real numbers $a, b, c$. Answer: $f(x)=\alpha x+\beta$ or $f(x)=\alpha x^{3}+\beta$ where $\alpha \in\{-1,0,1\}$ and $\beta \in \mathbb{R}$.
It is straightforward to check that the above functions satisfy the equation. Now let $f(x)$ satisfy the equation, which we denote $E(a, b, c)$. Then clearly $f(x)+C$ also does; therefore, we may suppose without loss of generality that $f(0)=0$. We start with proving Lemma. Either $f(x) \equiv 0$ or $f$ is injective. Proof. Denote by $\Theta \subseteq \mathbb{R}^{2}$ the set of points $(a, b)$ for which $f(a)=f(b)$. Let $\Theta^{*}=\{(x, y) \in \Theta$ : $x \neq y\}$. The idea is that if $(a, b) \in \Theta$, then by $E(a, b, x)$ we get 
$$ H_{a, b}(x):=\left(a b^{2}+b x^{2}+x a^{2}, a^{2} b+b^{2} x+x^{2} a\right) \in \Theta $$ 
for all real $x$. Reproducing this argument starting with $(a, b) \in \Theta^{*}$, we get more and more points in $\Theta$. There are many ways to fill in the details, we give below only one of them. Assume that $(a, b) \in \Theta^{*}$. Note that 
$$ g_{-}(x):=\left(a b^{2}+b x^{2}+x a^{2}\right)-\left(a^{2} b+b^{2} x+x^{2} a\right)=(a-b)(b-x)(x-a) $$ 
and 
$$ g_{+}(x):=\left(a b^{2}+b x^{2}+x a^{2}\right)+\left(a^{2} b+b^{2} x+x^{2} a\right)=\left(x^{2}+a b\right)(a+b)+x\left(a^{2}+b^{2}\right). $$ 
Hence, there exists $x$ for which both $g_{-}(x) \neq 0$ and $g_{+}(x) \neq 0$. This gives a point $(\alpha, \beta)=$ $H_{a, b}(x) \in \Theta^{*}$ for which $\alpha \neq-\beta$. Now compare $E(\alpha, 1,0)$ and $E(\beta, 1,0)$. The left-hand side expressions coincide, on the right-hand side we get $f(\alpha)-f\left(\alpha^{2}\right)=f(\beta)-f\left(\beta^{2}\right)$, respectively. Hence, $f\left(\alpha^{2}\right)=f\left(\beta^{2}\right)$ and we get a point $\left(\alpha_{1}, \beta_{1}\right):=\left(\alpha^{2}, \beta^{2}\right) \in \Theta^{*}$ with both coordinates $\alpha_{1}, \beta_{1}$ non-negative. Continuing squaring the coordinates, we get a point $(\gamma, \delta) \in \Theta^{*}$ for which $\delta>5 \gamma \geqslant 0$. Our nearest goal is to get a point $(0, r) \in \Theta^{*}$. If $\gamma=0$, this is already done. If $\gamma>0$, denote by $x$ a real root of the quadratic equation $\delta \gamma^{2}+\gamma x^{2}+x \delta^{2}=0$, which exists since the discriminant $\delta^{4}-4 \delta \gamma^{3}$ is positive. Also $x0$ for $x>0$ by injectivity. Assume that $f(x)$ for $x>0$ does not have the form $f(x)=x^{\tau}$ for a constant $\tau$. The known property of multiplicative functions yields that the graph of $f$ is dense on $(0, \infty)^{2}$. In particular, we may find positive $b1$. Also, such $b$ can be found if $f(x)=x^{\tau}$ for some $\tau0$ and some constant $\tau>0$. Dividing $E(a, b, c)$ by $(a-b)(b-$ $c)(c-a)=\left(a b^{2}+b c^{2}+c a^{2}\right)-\left(a^{2} b+b^{2} c+c^{2} a\right)$ and taking a limit when $a, b, c$ all go to 1 (the divided ratios tend to the corresponding derivatives, say, $\frac{a^{\tau}-b^{\tau}}{a-b} \rightarrow\left(x^{\tau}\right)_{x=1}^{\prime}=\tau$), we get $\tau^{3}=\tau \cdot 3^{\tau-1}, \tau^{2}=3^{\tau-1}, F(\tau):=3^{\tau / 2-1 / 2}-\tau=0$. Since function $F$ is strictly convex, it has at most two roots, and we get $\tau \in\{1,3\}$.