First, let me mention that log convexity of a function is implied by an analytic property, which appears to be more natural than log convexity itself. Namely, if $\mu$ is a Borel measure on $[0,\infty)$ such that the $r$th moment
$$f(r)=\int_{0}^{\infty}z^r d\mu(z)$$
is finite for all $r$ in the interval $I\subset \mathbb R$, then $\log f$ is convex on $I$.
Log convexity can be effectively used in derivation of various inequalities involving the gamma function (particularly, two-sided estimates of products of gamma functions). It is linked with the notion of Schur convexity which is itself used in many applications.
An appetizer. Let $m=\max x_i$, $s=\sum x_i$, $x_i > 0$, $i = 1,\dots,n$, then
$$[\Gamma(s/n)]^n\leq\prod\limits_{1}^{n}\Gamma (x_i)\leq \left[\Gamma\left(\frac{s-m}{n-1}\right)\right]^{n-1}\Gamma(m).\qquad\qquad\qquad (1)$$
(1) is trivial, of course, when all $x_i$ and $s/n$ are integers, but in general the bounds do not hold without assuming log convexity.
Edit added: a sketch of the proof. Let $f$ be a continuous positive function defined on an interval $I\subset \mathbb R$. One may show that the function $\phi(x)=\prod\limits_{i=1}^{n}f(x_i)$, $x\in I^n$ is Schur-convex on $I^n$ if and only if $\log f$ is convex on $I$. Thus the function
$$\phi(x)=\prod\limits_{i=1}^n \Gamma(x_i),\quad x_i>0,\qquad \quad\qquad\qquad\qquad\qquad\qquad\quad (2)$$
is Schur-convex on $I^n=(0,\infty)^n$. Since $x_i\le m$, $i=1,\dots,n$, and $\sum x_i=s$, it is easy to check that
$$x \prec \left(\frac{s-m}{n-1},\dots,\frac{s-m}{n-1},m\right).$$
The latter majorization and the fact that $\phi(x)$ defined by (2) is Schur-convex imply the upper bound (1). The lower bound follows from the standard majorization $x\succ (s/n,\dots,s/n)$.
Have a look at the recent short article by Marshall and Olkin concerning this and related inequalities for the gamma function.
Best Answer
I don't have a complete answer. As you say, many sources say that Euler did it, but Gronau gives compelling reason to doubt this. The best source I have found for this issue is "The early history of the factorial function" by Dutka, and for what it's worth I am convinced that Gronau's assessment is a fair one.
First, I'll summarize the usual story. Kline discusses this in chapter 19, section 5 of Mathematical Thought from Ancient to Modern Times (which falls in volume 2 of the paperback printing), and a more thorough source is Davis's article "Leonhard Euler's Integral: A Historical Profile of the Gamma Function". There is agreement in these sources that Euler solved the problem after unsuccessful attempts by Stirling, D. Bernoulli, and Goldbach, and that the first record of Euler's solution appears in outline form in a 1729 letter from Euler to Goldbach. This was expanded in subsequent letters and written up in the article to which Kristal Cantwell links (apparently the article was written in 1729 but not published until 1738). Euler's letters to Goldbach start on the third page of this pdf.
However, Gronau cites a letter by Bernoulli that was written a few days before Euler's and that contains at least a partial solution, possibly contradicting Kline and Davis. Dutka's paper goes into more detail and also claims that Euler's work was influenced by Bernoulli's earlier solution. I could only speculate on what led to the confusion among other authors, and I won't do so here. Perhaps it should be mentioned here (as is done by Gronau and Dutka) that Euler did much more than Bernoulli. For instance, Euler gave the first integral representations of the gamma function.
Edit: Because this answer is accepted and yet incomplete, I want to direct attention to Bruce Arnold's answer below. It contains a link to a copy of the too often neglected letter of D. Bernoulli cited by Gronau and Dutka.