To address your last comment on evaluating infinite sums of powers of logarithms:
Letting
$$(-1)^n \zeta^{(n)}(s)=\sum_{k=1}^\infty \frac{(\log k)^n}{k^s}$$
the question can be recast as the evaluation of the quantity $(-1)^n \zeta^{(n)}(0)$ for various $n\in \mathbb N$.
Apostol and Choudhury had given a bunch of formulae for evaluating $\zeta^{(n)}(s)$; in particular, there is the formula
$$(-1)^n \zeta^{(n)}(0)=\frac{\Im((-\log\,2\pi-i\pi/2)^{n+1})}{\pi(n+1)}+\frac1{\pi}\sum_{j=1}^{n-1} a_j j!\binom{n}{j}\Im((-\log\,2\pi-i\pi/2)^{n-j})$$
where
$$a_j=c_{j+1}+\sum_{\ell=0}^j \frac{(-1)^\ell \gamma_\ell}{\ell!} c_{j-\ell}$$
$$c_k=-\frac{\gamma c_{k-1}}{k}-\frac1{k}\sum_{\ell=1}^{k-1} (-1)^\ell \zeta(\ell+1) c_{k-\ell-1},\qquad c_0=1$$
the $\gamma_n$ are the Stieltjes constants, and $\gamma=\gamma_0$ is the Euler-Mascheroni Constant.
Here are a few explicit evaluations of the formula:
$$\zeta^{\prime\prime}(0)=\frac{\gamma^2}{2}-\frac{\pi^2}{24}-\frac{(\log (2\pi))^2}{2}+\gamma_1$$
$$-\zeta^{\prime\prime\prime}(0)=-\gamma^3-\frac32 \gamma ^2 \log(2\pi)+\frac{\pi^2}{8}\log(2\pi)+\frac{(\log(2\pi))^3}{2}-3\gamma\gamma _1-3\gamma_1\log (2\pi)-\frac32\gamma_2+\zeta(3)$$
As can be surmised, the closed forms become rather unwieldy as $n$ gets large...
An extension of these results to noninteger powers would involve having to appropriately define the differintegral of the zeta function; going the series route would now involve terms containing the incomplete gamma function, but I've no knowledge of any closed forms for the resulting sum.
A non rigorous approach, just simple "intuitive methode". Because the power is $x^2$, the period is 2. Make sure to include the zero point in it's polynomal notation. I will later try to write it for all power's of $x$.
$$ \sum_{x=1}^{\infty} (c)^{-x^2}$$
$$ \sum_{x=1}^{\infty} \sum_{n \in \mathbb{Z}}\ (\ln(c)x^2)^n(-1)^n/n!=\sum_{x=1}^{\infty} \sum_{n \in \mathbb{Z}}\ \frac{(x)^n (\ln(c))^{n/2}}{(n/2)!} (e^{i\pi n/2})\sum_{k=0}^1 \frac{e^{\frac{2i \pi*kn}{2}}}{2}$$
$$= \sum_{n \in \mathbb{Z}}\ \frac{\zeta(-n) (\ln(c))^{n/2}}{(n/2)!} \frac{(e^{i\pi n/2}+e^{3i\pi n/2})}{2}$$
First take the limit as n goes to -1. Which gives us the value $\frac{ \sqrt{\pi}}{2 \sqrt{\ln(c)}}$
There's another value at the "normal" order, mainly as $n$ goes to 0. At $n=0$ the value is $-\frac{1}{2}$.
Notice that both answers give a very very good approaximation for the sum, but are a very little off, (orders like $10^-{6}$, really small).
Mathematically I can't fully show why (yet). More on intuition based, let $h$ be a very small value, and you have rewrite the function you are trying to sum as $g(n)f(n)$ with $g(n)$ being a periodically zero, with an order $h$ near the zero so $g(n+h)$ and $f(n)$ being a function that grows more/equal then a constant value, the regularised (as in the upper or lower limit is infinity) the sum is not that well defined and you just have to be careful. Lots of horrible talk and explainations I soon will be ashamed of so I am hoping for a better explaination. So here's what I found to be the solution.
If $n$ is even, the zero's don't have any influence for the sum and can be regarded as 0.
If n is uneven this isn't the case:
$$= \sum_{n \in \mathbb{Z}}\ \frac{\zeta(-n) (\ln(c))^{n/2}}{(n/2)!} e^{i\pi n/2} \frac{(1+(-1)^n)}{2}$$
use that
$$\zeta(-n)=\frac{\zeta(n+1) n! (i^{n+1}+(-i)^{n+1})}{(2\pi)^{n+1}} $$
$$\frac{n!}{(n/2)!2^n}=\frac{((n-1)/2)!}{\sqrt{\pi}}$$
$$ \sum_{n \in \mathbb{Z}}\ \frac{\zeta(n+1)(\ln(c))^{n/2}\big(\frac{(n-1)}{2}\big)! }{\pi^{n+3/2}} \frac{(e^{i\pi n/2}+e^{3i\pi n/2}))((e^{i\pi (n+1)/2}+e^{3i\pi (n+1)/2}))}{4}=$$
We only had to look at when $n=2m-1$.
$$ \sum_{m \in \mathbb{Z}}\ \frac{\zeta(2m)(m-1)! (ln(c))^{m-1/2} }{2\pi^{2m+1/2}} i(e^{3i\pi m}-e^{i\pi (m)})(-1)^{m}$$
Because it's summed over all integers including the negative ones, we subsitute m=-m.
$$ \sum_{m \in \mathbb{Z}} \frac{\zeta(-2m)(-m-1)! \pi^{2m-1/2}}{2(\ln(c))^{m+1/2}} i(e^{-3i\pi m}-e^{-i\pi (m)}) (-1)^{m}$$
I've had the alternating/peridiocally zero's kinda horrible defined. But approach the limit of m=m+h at the zero's/infinities to look at the (uneven) zero's. So that $e^{-3i\pi (m+h)}-e^{-i\pi (m+h)}=-2hi\pi$ I know you should work this better out with it's period, and notice that everything is defined in the same "complex" period, maybe use trigonometry functions, but those functions are less intunitive. Again it's non rigorous, only to show that there's some way to a good solution with this methode.
And know that $(-n-1)!=h^{-1}/n!$ if n is an integer
$$ \sum_{m \in \mathbb{Z}} \frac{\zeta(-2m)(-m-1)! \pi^{2m+1/2}}{(\ln(c))^{m+1/2}} (h) (-1)^{m}$$
$$ \frac{ \sqrt{\pi}}{\sqrt{\ln(c)}}\sum_{m \in \mathbb{Z}} \frac{\zeta(-2m) \pi^{2m}}{m!(\ln(c))^{m}} (-1)^{m}$$
Officially I'd say you'd had to write it again out also as (and maybe you couldn't ignore the other zero's, I just don't know yet, it just works, I've to polish it obviously)
$$ \frac{ \sqrt{\pi}}{\sqrt{\ln(c)}}\sum_{m \in \mathbb{Z}} \frac{\zeta(-m) \pi^{m}}{(m/2)!(\ln(c))^{m/2}} (-1)^{m/2}(1+(-1)^m)$$
$$ \frac{ \sqrt{\pi}}{\sqrt{\ln(c)}}\sum_{x=1}^{\infty}\sum_{m \in \mathbb{Z}} \frac{x^{m} \pi^{m}}{(m/2)!(\ln(c))^{m/2}} (-1)^{m/2}(1+(-1)^m)$$
Which is the regularised sum of
$$\frac{ \sqrt{\pi}}{\sqrt{\ln(c)}}\sum_{x=1}^{\infty}e^{\frac{-\pi^2 x^2}{ln(c)}}$$
$$\sum_{x=1}^{\infty} (c)^{-x^2}=\frac{\sqrt{\pi}}{2 \sqrt{ln(c)}}-1/2+\frac{ \sqrt{\pi}}{\sqrt{\ln(c)}}\sum_{x=1}^{\infty}e^{\frac{-\pi^2 x^2}{\ln(c)}}.$$
The first two value's being your I(c^(-1),2) and Z(c^(-1),2) and the sum being your systematically error. And i can only agree that the systematically error is something counter intunitive at first, but also arise when you apply this summation methode to other series as i recall well known easier examples such as $\sum_{x=1}^{\infty} \frac{1}{2x^2-1}$.
In similair rough fashion for
$$ \sum_{x=1}^{\infty} (c)^{-x^d}$$
$$ \sum_{m \in \mathbb{Z}}\ \frac{(-1)^{1/d}\zeta(-m)(ln(c))^{m/d}}{(m/d)!} \sum_{k=0}^{d-1} \frac{e^{\frac{2ipi*km}{d}}}{d}$$
$$\sum_{m \in \mathbb{Z}}\ \frac{\zeta(-m)(ln(c))^{m/d}}{(m/d)!} \frac{e^{i\pi m/d}(e^{2 i \pi m}-1)}{d(e^{ 2i \pi m/d}-1)}$$
Now the value at m=0, will always be -1/2.
at limit m goes to -1 gives:
$$\frac{-\pi}{\sin(\pi/d) (\ln(c))^{1/d}(-\frac{1+d}{d})!}=\frac{(1/d)!}{\ln(c)^{1/d}}$$
Which is similair to the integeral you stated, if d is not even you got to indeed sum is over the other value's and don't forget the error.
A try of mine to find the error:
$$\sum_{m \in \mathbb{Z}} \frac{\zeta(-m)(ln(c))^{m/d}}{(m/d)!} \frac{e^{i\pi m/d}(e^{2 i \pi m}-1)}{d(e^{ 2i \pi m/d}-1)}$$
$$\sum_{m \in \mathbb{Z}/d\mathbb{Z}} \frac{-\zeta(-m)(ln(c))^{m/d}}{(m/d)!} \frac{h \pi}{d \sin(m \pi /d)}$$
$$\sum_{m \in \mathbb{Z}/d\mathbb{Z}} \frac{-\zeta(m+1) m!(ln(c))^{m/d}}{(2 \pi)^{m+1} (m/d)!} \frac{h \pi (i^{m+1}+(-i)^{m+1})}{d \sin(m \pi /d)}$$
$$\sum_{m \in \mathbb{Z}/d\mathbb{Z}} \frac{\zeta(m+1) m!(ln(c))^{m/d}}{(2 \pi)^{m} (m/d)!} \frac{h \pi (sin(m \pi/2)}{d \sin(m \pi /d)}$$
Use that:
$$\frac{m!}{(m/d)!}= d^{m+1/2} (2 \pi)^{(1-d)/2} \prod_{j=1}^{d-1} \big(\frac{m-j}{d}\big)!$$
$$\sum_{m \in \mathbb{Z}/d\mathbb{Z}} \frac{\zeta(m+1) (ln(c))^{m/d}}{} \frac{h (sin(m \pi/2)}{\sin(m \pi /d)}\frac{d^{m-1/2} (2 )^{(3-d)/2-m} (\pi)^{(1-d)/2-m}\prod_{j=1}^{d-1} \big(\frac{m-j}{d}\big)!}{}$$
$$\sum_{m \in \mathbb{Z}/d\mathbb{Z}} \frac{\zeta(-m+1) (ln(c))^{-m/d}}{} \frac{h (sin(-m \pi/2)}{\sin(-m \pi /d)}\frac{d^{-m-1/2} (2 )^{(3-d)/2+m} (\pi)^{(1-d)/2+m}\prod_{j=1}^{d-1} \big(\frac{-m-j}{d}\big)!}{}$$
Split up for the zero's in $\prod_{j=1}^{d-1} $
$$\sum_{m \in \mathbb{Z}/d\mathbb{Z}} \frac{\zeta(-m+1) (ln(c))^{-m/d}}{} \frac{h (sin(-m \pi/2)}{\sin(-m \pi /d)}\frac{d^{-m-1/2} (2 )^{(3-d)/2+m} (\pi)^{(3-d)/2+m}}{ \prod_{j=1}^{d-1}\sin(\frac{\pi (m+j)}{d}) \big(\frac{m+j}{d}-1\big)!}$$
$$\sum_{m \in \mathbb{Z}} \sum_{k=1}^{d-1}\frac{\zeta(-dm+1-k) (ln(c))^{(-dm-k)/d}}{} \frac{h (sin((-dm-k) \pi/2)}{\sin((-dm-k) \pi /d)}\frac{d^{(-dm-k)-1/2} (2 )^{(3-d)/2+dm+k} (\pi)^{(3-d)/2+dm+k}}{ \prod_{j=1}^{d-1}\sin(\frac{\pi (dm+j+k)}{d}) \big(\frac{ \pi (dm+j+k)}{d}-1\big)!}$$
$$\sum_{m \in \mathbb{Z}} \sum_{k=1}^{d-1}\frac{\zeta(-dm+1-k) (ln(c))^{(-dm-k)/d}}{m!} \frac{(sin((-dm-k) \pi/2)}{\sin((-dm-k) \pi /d)}\frac{d^{(-dm-k)-1/2} (2 )^{(3-d)/2+dm+k} (\pi)^{(1-d)/2+dm+k}}{ \prod_{(j=1)\wedge j\neq d-k}^{d-1}\sin(\frac{\pi (dm+j+k)}{d}) \big(\frac{dm+j+k}{d}-1\big)!}$$
$$\sum_{m \in \mathbb{Z}} \sum_{k=1}^{d-1}\frac{\zeta(-dm+1-k) (ln(c))^{(-dm-k)/d}}{m!} \frac{\sin((dm+k) \pi/2)}{\prod_{(j=0)\wedge j\neq d-k}^{d-1}\sin(\frac{\pi (dm+j+k)}{d})}\frac{d^{(-dm-k)-1/2} (2 )^{(3-d)/2+dm+k} (\pi)^{(1-d)/2+dm+k}}{ \prod_{(j=1)\wedge j\neq d-k}^{d-1} \big(\frac{dm+j+k}{d}-1\big)!}$$
$$\prod_{(j=0)\wedge j\neq d-k}^{d-1}\sin(\frac{\pi (dm+j+k)}{d})=d2^{1-d}(-1)^{?}$$
I don't have the minus sign for now. Will come back at that later, assuming it's not there cause i've already enough in one equation.
$$\sum_{m \in \mathbb{Z}} \sum_{k=1}^{d-1}\frac{\zeta(-dm+1-k) (ln(c))^{(-dm-k)/d}}{m!} \frac{\sin((dm+k) \pi/2)}{}\frac{d^{(-dm-k)-3/2} (2 )^{(1+d)/2+dm+k} (\pi)^{(1-d)/2+dm+k}}{ \prod_{(j=1)\wedge j\neq d-k}^{d-1} \big(\frac{dm+j+k}{d}-1\big)!}$$
to make it readable
$$\sum_{m \in \mathbb{Z}} \sum_{k=1}^{d-1}\frac{\zeta(-dm+1-k) (ln(c))^{(-dm-k)/d}}{m!} \frac{\sin((dm+k) \pi/2)}{\prod_{(j=1)\wedge j\neq d-k}^{d-1} \big(\frac{dm+j+k}{d}-1\big)!}\frac{(2 )^{(1+d)/2+dm+k} (\pi)^{(1-d)/2+dm+k}}{d^{(dm+k)+3/2} }$$
$$\sum_{m \in \mathbb{Z}} \sum_{k=1}^{d-1}\frac{\zeta(-dm+1-k) (ln(c))^{(-dm-k)/d}}{m!} \frac{\sin((dm+k) \pi/2)}{\prod_{(j=1)\wedge j\neq d-k}^{d-1} \big(\frac{dm+j+k}{d}-1\big)!} \frac{(\frac{2\pi}{d})^{dm+k+1/2} (\frac{2}{\pi})^{d/2}}{d}$$
$$\frac{\sqrt{\frac{2\pi}{d}} (\frac{2}{\pi})^{d/2}}{d} \sum_{m \in \mathbb{Z}} \sum_{k=1}^{d-1}\frac{\zeta(-dm+1-k)}{} \frac{ \bigg(\frac{2\pi}{d\sqrt[d]{ln(c)}}\bigg)^{(dm+k)}}{m!} \frac{\sin((dm+k) \pi/2)}{\prod_{(j=1)\wedge j\neq d-k}^{d-1} \big(\frac{dm+j+k}{d}-1\big)!} $$
I am already suprised the formula above seem to work for d=2.
If only the m didn't sneaked into the product life would be good.
Kinda used this is a note pad, i will clean it up one day, but results count.
$$\frac{\sqrt{\frac{2\pi}{d}} (\frac{2}{\pi})^{d/2}}{d} \sum_{m \in \mathbb{Z}} \sum_{k=1}^{d-1}\frac{\zeta(-dm+1-k)}{} \frac{ \bigg(\frac{2\pi}{d\sqrt[d]{ln(c)}}\bigg)^{(dm+k)}}{} \frac{\sin((dm+k) \pi/2)}{\prod_{(j=1)}^{d-1} \big(\frac{dm+j+k}{d}-1\big)!} $$
$$\frac{\sqrt{\frac{2\pi}{d}} (\frac{2}{\pi})^{d/2}}{d} \sum_{m \in \mathbb{Z}/d\mathbb{Z}} \frac{\zeta(-m+1)}{} \frac{ \bigg(\frac{2\pi}{d\sqrt[d]{ln(c)}}\bigg)^{(m)}}{} \frac{\sin((m) \pi/2)}{\prod_{(j=1)}^{d-1} \big(\frac{m+j}{d}-1\big)!} $$
It shows how I think you can proceed (I let it rest, as it got messy, and i probably made some errors, and made things unnessecary hard), with a doable solution of d=2, as was the original question. Maybe you can or should add $m \in d\mathbb{Z}$, to make it smooth, Again i do not know yet.
Best Answer
This is more of a (numerical) plausibility argument to support Gottfried's $\lim\limits_{n \to \infty} f(n) = \frac1{\log\;2}$ conjecture than an answer.
What struck me was that the error ratios $Q(n)$ as $n$ is doubled seem to indicate a decrease in error by a factor of 2. This reminded me of the decrease in error exhibited by the trapezoidal rule for numerical integration as the number of subintervals is doubled, leading me to think that Richardson extrapolation might be appropriate to use in numerically approximating the limit, just like the way Richardson extrapolation is used to improve trapezoidal estimates in Romberg quadrature.
Recall that what Richardson extrapolation essentially does is to fit an interpolating polynomial $p(x)$ to successive values of a sequence $s_i$, along with some auxiliary sequence $x_i$ that tends to 0 as $i\to\infty$ such that $p(x_i)=s_i$, and then take $p(0)$ to be the estimate of the limit of the sequence $s_i$. For this particular case, the error ratio appearing to approach 2 suggests that we take the auxiliary sequence $x_j=2^{-j}$.
The Richardson scheme proceeds as follows: letting $T_j^{(0)}=s_j$, we generate a triangular array of values
$$\begin{matrix}T_0^{(0)}&&&&\\T_1^{(0)}&T_1^{(1)}&&&\\T_2^{(0)}&T_2^{(1)}&T_2^{(2)}&&\\\vdots&\vdots&\vdots&\ddots&\\T_n^{(0)}&T_n^{(1)}&\cdots&\cdots&T_n^{(n)}\end{matrix}$$
whose entries are generated according to the recursive formula
$$T_j^{(n)}=T_{j}^{(n-1)}+\frac{T_{j}^{(n-1)}-T_{j-1}^{(n-1)}}{\frac{x_{j-n}}{x_j}-1}$$
and the $T_n^{(n)}$ (i.e., the "diagonal" elements of the array) form a sequence which (hopefully) converges faster to the limit than the original $s_j$. For this specific example, we have $s_j=f(2^{j+1})$ and $\frac{x_{j-n}}{x_j}=2^n$.
The Richardson scheme can be implemented such that only a one-dimensional array is needed, and this is what I will do in the Mathematica code that follows (specializing to the case $\frac{x_{j-n}}{x_j}=2^n$):
We then generate the (first fourteen members of the) sequence in the OP to 100 significant figures:
If we apply
richardson[]tovals, we get the following:and the result of the extrapolation matches $\frac1{\log\;2}$ to ~ 27 digits.
Generating more members of the sequence and computing with more digits might give a more accurate value, but instead of doing this, I decided to apply a different extrapolation algorithm, the Shanks transformation, with the reasoning that if two unrelated algorithms give (almost) the same answer for a given sequence, the answer is likely correct.
I will skip the details of the mathematics behind the Shanks transformation (but see this article) and instead present Mathematica code:
This is essentially the algorithm behind the hidden built-in Mathematica function
SequenceLimit[], but I have chosen to implement the algorithm explicitly for pedagogical reasons.We then have the results
with a result that matches $\frac1{\log\;2}$ to ~ 20 digits. (The slight loss of accuracy is the price for the generality of the Shanks transformation, which did not exploit the behavior of the error ratios.)
The results of these two different extrapolation algorithms lends credibility to the conjecture that the limit is $\frac1{\log\;2}$.