How to evaluate theta function's derivative numerically?

Numerical derivative

Based on $$f'(z_0)={1 \over 2\pi i}\,\int_\gamma {f(z) \, dz\over (z-z_0)^2}\,,$$ where $\gamma$ is a closed contour containing $z_0$ in its interior.

fPrime[z0_] :=
  1/8 Sum[f[z0 + dz]/dz,
          {dz, Exp[2 Pi I Most@Subdivide[0., 1., 8]]/1000}];

fPrime[0.1]
(*  -0.256724 + 1.47096 I  *)

---

Update:

Discretizing the integral with n = 2 points instead of n = 8 yields the central difference formula, and for a radius of Abs[dz] == 1*^-9, it will have a truncation error less than machine-precision for analytic functions whose higher-order derivatives do not grow too rapidly. To prevent round-off error overwhelming the truncation error, we compute f[z] at high precision. This is faster than the 8-point machine-precision code above on the OP's function (I suspect because SiegelTheta is somewhat expensive to compute). The 8-point formula with a radius of 1/1000 in fPrime has a relative error of $10^{-10}$ or less in a neighborhood of $z = 0.1 + 0i$. The function ND[] has a relative error of $10^{-5}$ or less. Over the square with ReIm[z] between ±1, the relative errors of fPrime and ND can be a couple of orders of magnitude larger, but fPrime2 below maintains machine-precision-accurate results.

ClearAll[f, fPrime2];
a = 1/10;
b = 2/10;
t = Exp[I 2 Pi/3];
f[z_] := SiegelTheta[{{a}, {b}}, {{t}}, z]
fPrime2[z_?NumericQ] := N@With[{z0 = SetPrecision[z, 32], r = 1*^-9},
    (f[r + z0] - f[-r + z0])/(2 r)
    ];

Symbolic derivative

For the OP's special case of SiegelTheta[], a symbolic derivative can be computed from the Sum[] of its theta series expansion, which returns a sum in terms of EllipticTheta[], whose derivative is implemented as EllipticThetaPrime[[]:

SiegelThetaPrime[{{a_}, {b_}}, {{t_}}, z_] = Simplify@D[
   Sum[Exp[
     I Pi ((n + {a}).{{t}}.(n + {a}) + 
        2 (n + {a}).(z + {b}))], {n, -Infinity, Infinity}],
   z]
(* 
(E^(-((I π (b + z)^2)/
  t)) π (-2 I (b + z) EllipticTheta[3, (π (b + a t + z))/t, 
     E^(-((I π)/t))] + 
   EllipticThetaPrime[3, (π (b + a t + z))/t, 
    E^(-((I π)/t))]))/(Sqrt[-I t] t)
*)

SiegelThetaPrime[{{1/10}, {1/5}}, {{Exp[I 2 π/3]}}, 0.1]
(*  -0.256724 + 1.47096 I  *)

You can compute a numerical derivative as follows

ClearAll[f, g];
Needs["NumericalCalculus`"]
a = 0.1;
b = 0.2;
t = Exp[I 2 Pi/3];
f[z_] := SiegelTheta[{{a}, {b}}, {{t}}, z]
g[z0_] := ND[f[z], z, z0]
g[0.1] 

(*-0.256725 + 1.47096 I*)

I haven't checked the result is correct