Double sum of reciprocal powers of integer LCMs   I am looking for approximations, or a closed form, if available, for  the sum  
$S(n,k)=\sum_{1\leq i,j,\leq n}  lcm(i,j)^{-k}$  
where lcm(i,j) is the least common multiple of integers $i,j$ and $k$ is a positive integer. 
For the case $k=1$, I can see that a very loose lower bound is
$3H_n-2$ where $H_n$ is the harmonic sum. Another lower bound is obtained by replacing $lcm(i,j)$ by
$ij$  which yields $H_n^2$ for $k=1$ and $H_n^{2k}$ in general. 
I am, however more interested in upper bounds, and small and even $k$ values 
I am certainly no expert in number theory, so my apologies if this is trivial, but at least the statement 
looks clean to me.
Remark: There seems to be some work on determinants of matrices with entries of the form above.
 A: For $k > 1$ the sums are bounded, by exactly $\zeta^3(k)/\zeta(2k),$ but maybe you're looking for something more exact than this? (The evaluation of these sums when $n\rightarrow\infty$ I thought was an old exercise in Polya and Szego's problem books, but apparently that's not the right source. I first found out about these identities from a friend, at any rate.) An outline of a proof goes as follows:
$$\zeta(2s)\sum_{i,j=1}^\infty \frac{1}{\textrm{lcm}(i,j)^s} = \zeta(2s)\sum_{(i,j)=1}\sum_{\lambda=1}^\infty \frac{1}{i^s j^s}\frac{1}{\lambda^s} = \sum_{u,v=1}^\infty\sum_{\lambda=1}^\infty \frac{1}{u^s v^s}\frac{1}{\lambda^s} = \zeta(s)^3$$
(Here $\lambda$ 'takes the role' of being the highest common factor of $i$ and $j$, in the first sum.)
I don't think there could be a very nice closed form for these sums, although you can collect terms and write them as a more standard Dirichlet series...
I'm curious, how does the question come up?
A: Thanks Brad.
Presumably, your ratio of zeta functions provide an upper bound as $n \rightarrow \infty$.
I had a quick look at Polya and Szego but did not manage to find it.
My actual sum is $n^k S(n,k)$ and the bound above would be $O(n^k)$ but I'd like a tighter upper bound, say $O(n^{k-1} \log(n))$ or something like this, if at all possible.
I am considering signal expansions along these kinds of regularly sampled support sets, and correlation moments between non-orthogonal expansion coefficients. 
If we let $n=5,$ and decompose $(f(1),\ldots,f(5))$ by taking inner products with the vectors 
$\{(1,1,1,1,1),(0,1,0,1,0),(0,0,1,0,0),\ldots,(0,\ldots,0,1)\}$ the correlations between these expansion coefficients lead directly to $n^k S(n,k).$
