We are currently looking for a fast, i.e. subquadratic, algorithm for the following equation:
$z\_m = \sum\_{i,j :\, (i \cdot j) = m} x\_i \cdot y\_j$.


That is, we are given two finite input vectors $x$ and $y$ and are trying to calculate 
the full vector $z$, where each element $z\_m$ in $z$ consists of the sum of the pairwise products of elements in 
the input vectors $x\_i$ and $y\_j$, whose indices $i$ and $j$ result in $m$ after 
multiplication.


A bit on the background: we are trying to efficiently multiply two probabilistic mass functions.


Since this problem is closely related to standard convolution (when changing product to 
sum):
$z\_m = \sum\_{i,j :\, (i + j) = m} x\_i \cdot y\_j$,


we already found out that it would be possible to log-transform the indices of the input 
vectors and rediscretize/bin them. Then you could use FFT-convolve to get the log-transformed
$z$ vector, which has to be exponentiated (on the indices) in the end. However, discretization most 
probably results in errors (depending on the resolution of the rediscretization) as well as increases the overall runtime because of the increasing length of the input vectors.


Does anyone know about an alternative solution to this problem (maybe with a hint to an 
existing implementation or on implementation details)? We have also looked at some 
theory on log-convolution/scale convolution (<https://en.wikipedia.org/wiki/Logarithmic_convolution>, with an implementation [here](http://casa.colorado.edu/~ajsh/FFTLog/)) since it seems to be at least very 
closely related. However, we are not sure if and how the rather complex theory applies to our discrete problem directly.



Subquadratic multiplication of probability mass functions (with log-convolution?)
There is little room for improvement over brute force because the average number of collisions is a very slowly growing function. 


The [Erdős multiplication table problem](http://oeis.org/A027424) asks how many numbers are products $ab$ with $1\le a,b \le n$. It is known that the count is $$\Theta\left(\frac{n^2}{(\log n)^c (\log\log n)^{3/2}}\right)$$ where $c=1-(1+\log\log 2)/\log 2 = 0.086$. So, if you multiply random variables with full support $\{1,...,n\}$, you get almost $n^2$ values with positive probability. The denominator is under $10$ for $n$ up to about $10^{18}$. 


