Motivation:
A more complicated generalization of Closed expressions for harmonic-like multiple sums are sums of a type for which the follwing is the simplest case
$$s=\sum _{m=1}^{\infty } \sum _{n=1}^{\infty } \frac{1}{m n \left(m^2+n^2\right)}\tag{1}$$
The difference to the previous type is the appearance of the sum of squares in the denominator.
Generating the denominator by
$$\frac{1}{m^2+n^2}=\int_0^1 x^{m^2+n^2-1} \, dx$$
the double sum factorizes under the integral and we have to consider this sum
$$s_1(x)=\sum _{k=1}^{\infty } \frac{x^{k^2}}{k}\tag{2}$$
for $0 \le x \le 1$. The numerical value at $x=1/2$ is
$$s_1(1/2) = 0.5319048623268223...$$
Now we have arrived at the
Question: Is there a closed form for $s_1(x)$, or at least for $s_1(1/2)$?
Possible first steps
If the denominator were 1 instead of $k$ the result would be $$s_2(x) =\sum _{k=1}^{\infty } x^{k^2}=\frac{1}{2} (\vartheta _3(0,x)-1)\tag{3}$$ which justifies the reference to the ellptic theta function in the question.
We could generate the other denominator by $\frac{1}{k}=\int_0^1 y^{k-1} \, dy$ and consider
$$s_3(x,y)=\sum _{k=1}^{\infty } x^{k^2} y^k$$
While writing this question I spotted in the automatic margin references this question Evaluate the double sum $\sum_{m=1}^{\infty}\sum_{n=1}^{m-1}\frac{ 1}{m n\left(m^2-n^2\right)^2}$ which looks similar. Maybe the techniques used there could lead to an answer.
We can reduce the motivating problem (1) to the calculation of an Euler sum with purely imaginary argument in the harmonic number.
Complex partial fraction of the denominator gives
$$\frac{1}{m^2+n^2}=\frac{1}{2 m}\left(\frac{1}{m-i n}+\frac{1}{m+i n}\right)$$
and the sum over $n$ can be done
$$\frac{1}{2} \sum _{n=1}^{\infty } \frac{1}{m^2 n (m-i n)}=\frac{H_{i m}}{2 m^3}$$
Hence we are left with the Euler sum with imaginary argument
$$s_4 = \sum _{m=1}^{\infty } \frac{H_{i m}}{m^3}\tag{4}$$
where as usual $H_z$ is the harmonic number of $z$, and we have
$$s = \Re s_4\tag{5}$$
