In this discussion a program for computing first terms of the series expansion for $\sigma$-function associated with a genus two curve is presented. The series is computed for the curve in its canonical form. $\sigma$-Functions associated with curves in other forms are obtained by a proper transformation of parameters of the curve, and coordinates of the Jacobian varieties.
Canonical curve
Let the canonical form of genus two curves be defined by the equation $$\mathcal{C}_{(2,5)}:\quad f(x,y;\lambda) \equiv - y^2 + x^5 + \lambda_4 x^3 + \lambda_6 x^2 + \lambda_8 x + \lambda_{10} = 0. \tag{1} $$ The canonical form is the $(2,5)$-curve defined in [1], which is the Weierstass canonical form described in [2], $\S\S$ 60–63. Let $\lambda \equiv (\lambda_4, \lambda_6, \lambda_8, \lambda_{10} ) $ be parameters of the canonical curve (1). Let first and second kind differentials be defined by
$$\begin{aligned} &\mathrm{d} u= \begin{pmatrix} \mathrm{d} u_1 \\ \mathrm{d} u_3 \end{pmatrix} = \begin{pmatrix} x \\ 1 \end{pmatrix} \frac{\mathrm{d}x}{\partial_y f(x,y;\lambda)},\\ &\mathrm{d} r = \begin{pmatrix} \mathrm{d} r_1 \\ \mathrm{d} r_3 \end{pmatrix} = \begin{pmatrix} x^2 \\ 3 x^3 + \lambda_4 x \end{pmatrix} \frac{\mathrm{d}x}{\partial_y f(x,y;\lambda)}. \end{aligned} $$
The $\sigma$-function is obtained as a solution of the system of heat equations
$$Q_k \sigma(u_1,u_3;\lambda) = 0 $$ with the initial condition given by the Weierstrass--Schur polynomial $u_3 - \frac{1}{3} u_1^3$.
Operators $Q_k$, which annihilate the $\sigma$-function, have the following form, see arXiv:1509.01490, Sect. 4 p.6, and also arXiv:1711.08395, Sect 6.4,
$$\begin{aligned} &Q_0 = - u_1 \partial_{u_1} − 3 u_3 \partial_{u_3} + 3 + L_0,\\ &Q_2 = − \tfrac{1}{2} \partial_{u_1,u_1} + \tfrac{4}{5} \lambda_4 u_3 \partial_{u_1} − u_1 \partial_{u_3} + \tfrac{3}{10} \lambda_4 u_1^2 − \tfrac{1}{10}(15 \lambda_8 − 4 \lambda_4^2) u_3^2 + L_2,\\ &Q_4 = − \partial_{u_1,u_3} + \tfrac{6}{5} \lambda_6 u_3 \partial_{u_1} − \lambda_4 u_3 \partial_{u_3} + \tfrac{1}{5} \lambda_6 u_1^2 − \lambda_8 u_1 u_3 − \tfrac{1}{10} (30 \lambda_{10} − 6 \lambda_6 \lambda_4) u_3^2 + \lambda_4 + L_4,\\ &Q_6 = − \frac{1}{2} \partial_{u_3,u_3} + \tfrac{3}{5} \lambda_8 u_3 \partial_{u_1} + \tfrac{1}{10} \lambda_8 u_1^2 − 2 \lambda_{10} u_1 u_3 + \tfrac{3}{10} \lambda_8 \lambda_4 u_3^2 + \tfrac{1}{2} \lambda_6 + L_6, \end{aligned} $$ where the linear differential operators with respect to parameters $\lambda$ are obtained by
$$ \begin{pmatrix} L_0 \\ L_2 \\ L_4 \\ L_6 \end{pmatrix} = \begin{pmatrix} 4 \lambda_4 & 6 \lambda_6 & 8 \lambda_8 & 10 \lambda_{10} \\ 6 \lambda_6 & 8 \lambda_8 − \tfrac{12}{5} \lambda_4^2 & 10 \lambda_{10} − \tfrac{8}{5} \lambda_6 \lambda_4 & − \tfrac{4}{5} \lambda_8 \lambda_4 \\ 8 \lambda_8 &10 \lambda_{10} − \tfrac{8}{5} \lambda_6 \lambda_4 & 4\lambda_8 \lambda_4 − \tfrac{12}{5} \lambda_6^2 & 6 \lambda_{10} \lambda_4 − \tfrac{6}{5} \lambda_8 \lambda_6 \\ 10 \lambda_{10} & − \tfrac{4}{5} \lambda_8 \lambda_4 & 6 \lambda_{10} \lambda_4 − \tfrac{6}{5} \lambda_8 \lambda_6 & 4\lambda_{10} \lambda_6 − \tfrac{8}{5} \lambda_8^2 \end{pmatrix} \begin{pmatrix} \partial_{\lambda_4} \\ \partial_{\lambda_6} \\ \partial_{\lambda_8} \\ \partial_{\lambda_{10}} \end{pmatrix}. $$ The system is solved by means of representing the unknown $\sigma$-function as a series with terms formed according to the Sato weight, see arXiv:1711.08395, Sect 6.4. Recurrence relations and the corresponding series expansion are implemented below.
As an example, $\sigma$-function is computed with different accuracy: up to the weights $-13$, and $-25$ with respect to coordinates $u_1$, $u_3$ of the Jacobian variety of $\mathcal{C}_{(2,5)}$. Fundamental cubic relations and the identities for 4-index $\wp$-functions are used for verification.
A series representation is computed only for the $\sigma$-function associated with the canonical curve $\mathcal{C}_{(2,5)}$. This is enough to obtain the $\sigma$-function associated with any curve of genus two. A transformation of the $\sigma$-function along with the corresponding transformation between curves is illustrated below by two examples: the case of a more general $(2,5)$-curve with term $\mu_2 x^4$, and the case of a $(2,6)$-curve.
(2,5)-Curve with additional term $\mu_2 x^4$
Let a genus two curve be defined by the equation $$\tilde{\mathcal{C}}_{(2,5)}:\quad \tilde{f}(\tilde{x},\tilde{y};\mu) \equiv - \tilde{y}^2 + \tilde{x}^5 + \mu_2 \tilde{x}^4+ \mu_4 \tilde{x}^3 + \mu_6 \tilde{x}^2 + \mu_8 \tilde{x} + \mu_{10} = 0. \tag{2} $$ The curve (1) is obtained from (2) by the transformation $\tilde{x} = x - \tfrac{1}{5} \mu_2$, $\tilde{y} = y$. The backward transformation $$x = \tilde{x} + \tfrac{1}{5} \mu_2, \qquad y = \tilde{y}. \tag{3} $$ along with the transformation of parameters $\lambda$ into $\mu$
$$\begin{aligned} &\lambda_{10} = \tfrac{4}{5^5}\mu _2^5 - \tfrac{1}{5^3} \mu _4 \mu_2^3 +\tfrac{1}{5^2} \mu_6 \mu_2^2-\tfrac{1}{5} \mu_8 \mu_2+\mu_{10}\\ &\lambda_8 = -\tfrac{3}{5^3}\mu_2^4 + \tfrac{3}{5^2} \mu_4 \mu_2^2 - \tfrac{2}{5}\mu_6 \mu_2 + \mu_8\\ &\lambda_6 = \tfrac{4}{5^2} \mu_2^3 -\tfrac{3}{5}\mu_4 \mu_2 + \mu_6\\ &\lambda_4 = \mu_4-\tfrac{2}{5} \mu_2^2, \end{aligned} \tag{4}. $$ induces the transformation of first and second kind differentials $$\begin{pmatrix} \mathrm{d} u \\ \mathrm{d} r \end{pmatrix} = \begin{pmatrix} A & B \\ C & D \end{pmatrix} \begin{pmatrix} \mathrm{d} \tilde{u} \\ \mathrm{d} \tilde{r} \end{pmatrix}. \tag{5} $$ where $$\begin{aligned} &A = \begin{pmatrix} 1 & \tfrac{1}{5} \mu_2 \\ 0 & 1 \end{pmatrix},& &B = 0,& \\ &C = \begin{pmatrix} \tfrac{2}{5} \mu_2 & \tfrac{1}{25} \mu_2^2 \\ - \tfrac{1}{25} \mu_2^2 & \tfrac{1}{5} \mu_2 (\mu_4 - \tfrac{7}{25} \mu_2^2) \end{pmatrix},& &D = \begin{pmatrix} 1 & 0 \\ -\tfrac{1}{5} \mu_2 & 1 \end{pmatrix}.& \end{aligned} $$
Then the $\sigma$-function associated with $\tilde{\mathcal{C}}_{(2,5)}$ is computed by $$ \tilde{\sigma} (\tilde{u}; \mu) = \exp \Big( \tfrac{1}{2} \tilde{u}^t T \tilde{u} \Big) \sigma \big(A \tilde{u}; \lambda(\mu)\big), $$ where $T$ is a symmetric matrix, $$T = D^{-1} C = \begin{pmatrix} \tfrac{2}{5} \mu_2 & \tfrac{1}{25} \mu_2^2 \\ \tfrac{1}{25} \mu_2^2 & \tfrac{1}{5} \mu_2 (\mu_4 - \tfrac{6}{25} \mu_2^2) \end{pmatrix}, $$ and $\lambda(\mu)$ is given by (4). This series coincides with the one given in arXiv:1711.08395, p.48, and also in [3], p.89.
(2,6)-Curve
Now, we consider a genus two curve defined by the equation $$\begin{split} \bar{\mathcal{C}}_{(2,6)}:\quad \bar{f}(\bar{x},\bar{y};\mu) \equiv & - \bar{y}^2 + P(\bar{x}) \\ &- \bar{y}^2 + \nu_0 \bar{x}^6 + \nu_1 \bar{x}^5 + \nu_2 \bar{x}^4 + \nu_3 \bar{x}^3 + \nu_4 \bar{x}^2 + \nu_5 \bar{x} + \nu_6 = 0. \end{split} \tag{6} $$ Let $e_0$, such that $P(e_0)=0$, maps to infinity on $\tilde{\mathcal{C}}_{(2,5)}$. Then the transformation from $\bar{\mathcal{C}}_{(2,6)}$ to $\tilde{\mathcal{C}}_{(2,5)}$ has the form $$\bar{x} = e_0 + \frac{c}{\tilde{x}},\quad \bar{y} = \frac{c^3 \tilde{y}}{h \tilde{x}^3},\quad \text{where}\quad h^2 = \frac{c^5}{P'(e_0)}, $$ and $c$ is an arbitrary constant. This transformation is adopted from [4], p.144.
The backward transformation from $\tilde{\mathcal{C}}_{(2,5)}$ to $\bar{\mathcal{C}}_{(2,6)}$ $$\tilde{x} = \frac{c}{\bar{x}-e_0},\quad \tilde{y}=\frac{h \bar{y}}{(\bar{x}-e_0)^3}, \quad \mu_{2k} = \frac{c^k P^{k+1}(e_0)}{(k+1)! P'(e_0)} \tag{7} $$ induces the transformation $$\begin{pmatrix} \mathrm{d} \tilde{u} \\ \mathrm{d} \tilde{r} \end{pmatrix} = \begin{pmatrix} \tilde{A} & \tilde{B} \\ \tilde{C} & \tilde{D} \end{pmatrix} \begin{pmatrix} \mathrm{d} \bar{u} \\ \mathrm{d} \bar{r} \end{pmatrix}. \tag{8} $$ where $$\begin{aligned} &\tilde{A} = -\frac{c}{h} \begin{pmatrix} 0 & c \\ 1 & -e_0 \end{pmatrix},& &\tilde{B} = 0,& \\ &\tilde{C} = \frac{h}{c^2} \begin{pmatrix} 2 \nu_0 e_0^3 + \nu_1 e_0^2 & 4 \nu_0 e_0^4 + 3 \nu_1 e_0^3 + 2 \nu_2 e_0^2 + \nu_3 e_0 \\ c (6 \nu_0 e^2 + 2 \nu_1 e_0) & - c (4 \nu_0 e_0^3 + \nu_1 e_0^2) \end{pmatrix},& &\tilde{D} = - \frac{h}{c^2} \begin{pmatrix} e_0 & 1 \\ c & 0 \end{pmatrix}.& \end{aligned} $$
The $\sigma$-function associated with $\bar{\mathcal{C}}_{(2,6)}$ is computed by $$\bar{\sigma} (\bar{u}; \nu ) = \exp \Big( \tfrac{1}{2} \bar{u}^t \tilde{T} \bar{u} \Big) \tilde{\sigma} \big(\tilde{A} \bar{u}; \mu(\nu)\big), $$ where $\tilde{T}$ is a symmetric matrix, $$\tilde{T} = \tilde{D}^{-1} \tilde{C} = \begin{pmatrix} -6\nu_0 e_0^2 - 2 \nu_1 e_0 & 4 \nu_0 e_0^3 + \nu_1 e_0^2 \\ 4 \nu_0 e_0^3 + \nu_1 e_0^2 & -8\nu_0 e_0^4 - 4 \nu_1 e_0^3 - 2 \nu_2 e_0^2 - \nu_3 e_0 \end{pmatrix}, $$ and $\mu(\nu)$ is given by (7).
In order to find the expression for $\bar{\sigma}$ through $\sigma$, we combine the transformation from $\tilde{\mathcal{C}}_{(2,5)}$ to $\bar{\mathcal{C}}_{(2,6)}$ with transformation from $\tilde{\mathcal{C}}_{(2,5)}$ to $\mathcal{C}_{(2,5)}$, and obtain $$x = \frac{c}{\bar{x}-e_0} + \frac{1}{5} \frac{c P''(e_0)}{2! P'(e_0)},\quad y = \frac{h \bar{y}}{(\bar{x}-e_0)^3}. $$ Then, $\lambda(\nu)$ is obtained by combining (4) with $\mu(\nu)$ given by (7). Next,
$$ \begin{pmatrix} \mathrm{d} u \\ \mathrm{d} r \end{pmatrix} = \begin{pmatrix} A \tilde{A} + B \tilde{C} & A \tilde{B} + B \tilde{D} \\ C \tilde{A} + D \tilde{C} & C \tilde{B} + D \tilde{D} \end{pmatrix} \begin{pmatrix} \mathrm{d} \bar{u} \\ \mathrm{d} \bar{r} \end{pmatrix}. \tag{9} $$ and so, taking into account that $B=0$, and $\tilde{B}=0$, we have $$\bar{\sigma} (\bar{u}; \nu ) = \exp \Big( \tfrac{1}{2} \bar{u}^t \bar{T} \bar{u} \Big) \sigma \big(A \tilde{A} \bar{u}; \lambda(\nu)\big), $$ where $\bar{T} = (D \tilde{D})^{-1} (C \tilde{A} + D \tilde{C} ) = \tilde{D}^{-1} T \tilde{A} + \tilde{T}. $
References
[1] Buchstaber,V. M., Enolskii,V.Z., Leikin,D.V., Rational analogs of abelian functions, Functional Analysis and Its Applications, 33:2 (1999) 83--94.
[2] Baker, H.F., Abel's theorem and the allied theory of theta functions, Cambridge Univ. Press, Cambridge, 1897.
[3] Baker, H.F., Multiply periodic functions, Cambridge Univ. Press, Cambridge, 1907.
[4] Baker, H.F., On a system of differential equations leading to periodic functions, Acta Mathematica, 27 (1903) pp. 135--156.