How find all $f:\mathbb R\to\mathbb R$ such that $f\bigl(x\cdot f(y)\bigr)=y\cdot f(x)+kxy$

Theorem. If $k<-\frac14$, no such function exists. If $k=-\frac14$, the only solution is $f(x)=\frac12x$. If $k> -\frac14$ and $k\ne 0$ there are exactly two solutions $f(x)=cx$ with $c=\frac{1\pm\sqrt{1+4k}}{2}$. If $k=0$ there are exactly two continuous solutions $f(x)=0$ and $f(x)=x$ and infinitely many non-continuous solutions (explicitly described in the proof below and most of them obtainable only using the Axiom of Choice).

Proof: The function $f$ might be linear, $f(x)=cx$. In that case we obtain $c^2xy=cxy+kxy$ for all $x,y$, i.e. $c=\frac{1\pm\sqrt{1+4k}}{2}$ (only possible if $k\ge -\frac14$, of course). One checks that these are indeed solutions. So from now on assume that $f$ is not linear.

As $f$ is not linear, there exixts $x_0$ with $f(x_0)\ne -kx_0$. Let $y_0=x_0f\left(\frac1{f(x_0)+kx_0}\right)$. Then $$f(y_0) = f\Biggl(x_0f\left(\frac1{f(x_0)+kx_0}\right)\Biggr)=\frac1{f(x_0)+kx_0}\bigl(f(x_0)+kx_0\bigr)=1.$$ Then $$\tag1\label1 f(x)=y_0f(x)+ky_0 x\quad\text{for all }x.$$ If $y_0\ne1$, we can solve for $f(x)$ and obtain that $f(x)=\frac{ky_0}{1-y_0}\cdot x$ is linear, contrary to assumption. Therefore $y_0=1$ and hence (with plugging in $x=1$ into \eqref{1}) we find $k=0$.

Starting all over again with $k=0$, the functional equation becomes $$\tag2\label2 f\bigl(xf(y)\bigr)=yf(x)\quad\text{for all }x,y$$ and has a more or less obvious solution $$\tag3\label3 f(x)=\begin{cases}\tfrac1x&\text{if } x\ne 0\\0&\text{if }x=0.\end{cases}$$ More generally, let $\mathbb R=U\oplus V$ be any direct sum decomposition of $\mathbb R$ as $\mathbb Q$ vector space and let $g\colon \mathbb R\to\mathbb R$ be given by $g(u+v)=u-v$ for $u\in U, v\in V$. Then $$\tag4\label4 f(x)=\begin{cases}0&\text{if }x=0,\\ e^{g(\ln x)}&\text{if }x>0, \\-f(-x)&\text{if }x<0\end{cases}$$ is a solution (check!) that is so to speak combined from $f(x)=x$ and $f(x)=\frac1x$. Note that without invoking the Axiom of Choice we know only two such direct sum decompositions, namely those with one summand zero. These lead to the identity function ad \eqref{3}, respectively.

Can one show that all solutions for $k=0$ are of this form \eqref{4}? Not quite.

As $f$ is not linear, there exists $y$ with $f(y)\ne 0$. Then for all $t\in \mathbb R$ there exists $x$ with $t=xf(y)$. We find that $f\bigl(f(t)\bigr)=f\Bigl(f\bigl(xf(y)\bigr)\Bigr)=f\bigl(yf(x)\bigr)=xf(y)=t$, i.e. $f$ is an involution and especially $f$ is bijective. Then substituting $f(y)$ for $y$ in \eqref{2}, we obtain $$f(xy)=f(x)f(y)$$ and especially $f\restriction_{\mathbb R^\times}$ is a group automorphism. We have $f(1)=1$ and $f(-1)=-1$ (because $f(-1)^2=f(1)$ and $f(-1)\ne f(1)$) and therefore $f(-x)=-f(x)$ for all $x$. Let $$\hat f(x)=\operatorname{sgn}(x)|f(x)|.$$ Then one checks that $ \hat f(xy)=\hat f(x)\hat f(y)$ and $\hat f\left(\hat f(x)\right)=x$ and hence $\hat f$ is also a solution of the original functional equation, but with the additional property that $\hat f(x)>0$ iff $x>0$. Given $x>0$, we have $x\hat f(x)>0$ and can let $u=\sqrt{x\hat f(x)}$, $v=\frac xu$ so that $x=uv$. One checks that $\hat f\left(u^2\right)=u^2$ and hence $\hat f(u)=u$. Also, $\hat f(v)=\frac{\hat f(x)}{\hat f(u)}=\frac 1v$. After taking logarithms, this corresponds exactly to the direct sum decomposition used above to construct \eqref{4}. We conclude that (for $k=0$) for any nonzero $f$ at least the sign-corrected function $\hat f$ is of the form \eqref{4}.

Remains to check how strange the signs of $f$ can be. Since $f$ and $\hat f$ are both multiplicative, so is their quotient on $\mathbb R^\times$, which induces a group homomorphism from $\mathbb R^\times/\{\pm1\}\to\{\pm1\}$. And conversely, any $\hat f$ as given by \eqref{4} for a suitable direct sum decomposition of $\mathbb R$, together with a homomorphism $\mathbb R^\times/\{\pm1\}\to\{\pm1\}$ (essentially given by a subgroup of index $2$ in $\mathbb R_{>0}$ in the nontrivial case; again, without the Axiom of Choice we know only the trivial homomorphism) gives us a solution $f$ for the original functional equation. $_\square$