ConformalLabpp/doc/math/hyperideal-hessian-derivation.md

Analytic Hessian of the Hyper-Ideal Discrete Conformal Energy

Status: research note, mathematical preparation for Phase 9b-analytic (see doc/roadmap/research-track.md). The current shipped Hessian (hyper_ideal_hessian_block_fd, v0.9.0) is block finite-difference; this note derives the closed-form replacement via the Schläfli identity and the chain rule through the building blocks zeta, zeta13, zeta14, zeta15, lij, alpha_ij, sigma_i, sigma_ij declared in code/include/hyper_ideal_geometry.hpp and assembled by face_angles_from_local_dofs(...) in code/include/hyper_ideal_functional.hpp.

Target reader. A mathematician fluent with Springborn (2020) and the Bobenko-school discrete-conformal apparatus; the goal is verifiability of each derivative line against the source papers.

Conventions.

• b_i \in \mathbb{R} — log scale factor at hyper-ideal vertex i (the DOF). An "ideal" vertex has no b_i DOF; we encode that by vi_var = false.
• a_{ij} \in \mathbb{R} — intersection-angle DOF on edge ij (in the Springborn $a$-parametrisation, edges may be conventional, ideal–ideal, or hyper-ideal–hyper-ideal).
• \ell_{ij} — effective hyperbolic length of edge ij in the auxiliary truncated tetrahedron.
• \beta_i^{(f)} — interior angle of the auxiliary hyperbolic triangle on face f at vertex i.
• \alpha_{ij}^{(f)} — dihedral angle of the truncated tetrahedron at edge ij, contributed by face f.
• \Theta_v, \theta_e — prescribed cone / intersection-angle targets.
• V(f) — hyperbolic volume of the truncated tetrahedron associated with face f.

---

1. The Schläfli identity and the gradient

1.1 Energy

The Springborn (2020, §4) hyper-ideal energy on a triangulated surface M with DOF vector x = (b, a) is

E(b, a) \;=\; \sum_{f \in F} U(f) \;-\; \sum_{e \in E} \theta_e\, a_e \;-\;
                                     \sum_{v \in V} \Theta_v\, b_v ,

where the per-face contribution is

U(f) \;=\; \sum_{e \in f} a_e\, \alpha_e^{(f)} \;+\;
           \sum_{i \in f} b_i\, \beta_i^{(f)} \;+\; 2\, V(f) .

(Each face contributes its three edges and three vertices; ideal vertices contribute b_i = 0 trivially.) This matches the face_energy(...) routine in hyper_ideal_functional.hpp line by line.

1.2 First-order Schläfli (gradient)

The classical Schläfli differential identity (Schläfli 1858/60; Milnor 1982; Vinberg 1993, Ch. 7) for the volume of any compact hyperbolic polyhedron P \subset \mathbb{H}^3 reads

\boxed{\;\; -2\, dV \;=\; \sum_{e \subset P} \ell_e \, d\alpha_e \;\;}

in the closed (finite-vertex) case, where \ell_e, \alpha_e are edge length and dihedral angle. For the truncated / hyper-ideal extension we follow Bobenko–Springborn–Schief and write the identity in the mixed form that respects the truncation:

\boxed{\;\; 2\, dV \;=\; \sum_e a_e\, d\alpha_e \;+\; \sum_v b_v\, d\beta_v
                       \;-\; \sum_e \alpha_e\, da_e \;-\; \sum_v \beta_v\, db_v . \;\;}
\tag{S1}

Combined with the trivial d(a_e \alpha_e) = a_e\, d\alpha_e + \alpha_e\, da_e identity, (S1) is equivalent to

d\!\left(\sum_e a_e \alpha_e + \sum_v b_v \beta_v + 2V\right)
   \;=\; \sum_e \alpha_e\, da_e \;+\; \sum_v \beta_v\, db_v ,

i.e. dU(f) = \sum \alpha_e\, da_e + \sum \beta_v\, db_v on each face. Summing over faces and subtracting the linear \theta, \Theta terms,

\frac{\partial E}{\partial b_v} \;=\; \Big(\!\!\sum_{f \ni v} \beta_v^{(f)}\Big) - \Theta_v,
\qquad
\frac{\partial E}{\partial a_e} \;=\; \Big(\!\!\sum_{f \supset e} \alpha_e^{(f)}\Big) - \theta_e ,

which is precisely the gradient implemented in evaluate_hyper_ideal(...) (pass 3).

1.3 Second-order Schläfli (Hessian)

Differentiating (S1) once more (Cho–Kim 1999, Lemma 2.1; Rivin 1994) yields the second Schläfli identity

\boxed{\;\; 0 \;=\; \sum_e da_e \wedge d\alpha_e \;+\; \sum_v db_v \wedge d\beta_v . \;\;}
\tag{S2}

Since the wedge product is antisymmetric in the differential factors, (S2) forces the bilinear form

H(f) := \begin{pmatrix}
   \partial \beta_i / \partial b_j  & \partial \beta_i / \partial a_e \\
   \partial \alpha_e / \partial b_j & \partial \alpha_e / \partial a_{e'}
\end{pmatrix}

to be symmetric on each face (and so on the whole mesh). This is the deep reason that one may compute the Hessian

H \;=\; \frac{\partial^2 E}{\partial x^2} \;=\;
\frac{\partial (\beta - \Theta,\; \alpha - \theta)}{\partial (b, a)}

by accumulating only the upper triangle face by face, and that the mesh-level matrix inherits PSD-ness from the per-face 6 \times 6 blocks (Springborn 2020 §4.3).

1.4 Strategy

The full chain through face_angles_from_local_dofs(...) is

(b_i, a_e) \;\xrightarrow{\;\text{lij}\;}\; \ell_e
            \;\xrightarrow{\;\zeta\;}\; \beta_i
            \;\xrightarrow{\;\zeta,\sigma\;}\; \alpha_e .

We derive each arrow in Sections 2–3 and assemble the 6 \times 6 face Jacobian in Sections 4–5. Symmetry is then a consequence of (S2) and a useful numerical check.

---

2. Derivative of \zeta (interior angle from edge lengths)

2.1 Definition

For an auxiliary hyperbolic triangle with side lengths x, y, z (opposite to vertices X, Y, Z),

\zeta(x, y, z) \;=\; \arccos\!\left(\frac{\cosh x \cosh y - \cosh z}{\sinh x \sinh y}\right) ,

returning the interior angle at the vertex between the sides of length x and y (so z is opposite). This is the hyperbolic law of cosines solved for the included angle, and it is exactly zeta(...) in hyper_ideal_geometry.hpp.

2.2 Partial derivatives

Let

N := \cosh x \cosh y - \cosh z,
\qquad D := \sinh x \sinh y,
\qquad N / D = \cos \beta .

Then \sin \beta \cdot d\beta = -d(N/D) = (D\, dN - N\, dD)/D^2 \cdot (-1). Computing the partials,

\partial_x N = \sinh x \cosh y,\qquad
\partial_y N = \cosh x \sinh y,\qquad
\partial_z N = -\sinh z ,

\partial_x D = \cosh x \sinh y,\qquad
\partial_y D = \sinh x \cosh y,\qquad
\partial_z D = 0 .

Hence

\partial_x (\cos \beta) \;=\; \frac{(\sinh x \cosh y)\,\sinh x \sinh y - (\cosh x \cosh y - \cosh z)\, \cosh x \sinh y}{\sinh^2 x \sinh^2 y}

= \;\frac{\sinh y \big[ \sinh^2 x \cosh y - \cosh x (\cosh x \cosh y - \cosh z) \big]}{\sinh^2 x \sinh^2 y}
\;=\;\frac{\cosh x \cosh z - \cosh y}{\sinh^2 x \sinh y}

(using \sinh^2 x - \cosh^2 x = -1, which gives the cancellation \sinh^2 x \cosh y - \cosh^2 x \cosh y = -\cosh y).

Combined with -\sin \beta\, \partial_x \beta = \partial_x (\cos\beta) this yields the hyperbolic dual law of cosines in derivative form (Cho–Kim 1999, Lemma 3.1):

\boxed{\;\; \frac{\partial \beta}{\partial x}
   \;=\; \frac{\cosh y - \cosh x \cosh z}{\sin \beta \cdot \sinh^2 x \sinh y}
   \;=\; -\frac{\cos \beta_y}{\sinh x \sin \beta} \,\cdot\, \frac{1}{?} \;\;}

We prefer the normalised form (sine rule). Recall the hyperbolic sine rule on a triangle \beta, \beta_y, \beta_z: \sinh x / \sin \beta_x' = \sinh y / \sin \beta_y' = \sinh z / \sin \beta_z', where \beta_x' is opposite to side x. In our convention \beta = \beta_z' (the angle opposite z — wait: \zeta(x,y,z) returns the angle between sides x and y, i.e. opposite to side $z$), so let \beta := \beta_z', \beta_x' := angle opposite x, \beta_y' := angle opposite y. Then by the dual law of cosines applied to side y, \cosh y = \cosh x \cosh z - \sinh x \sinh z \cos \beta_y', whence \cosh y - \cosh x \cosh z = -\sinh x \sinh z \cos \beta_y'. Substituting,

\frac{\partial \beta}{\partial x}
   \;=\; \frac{-\sinh x \sinh z \cos \beta_y'}{\sin\beta \cdot \sinh^2 x \sinh y}
   \;=\; \frac{-\sinh z \cos \beta_y'}{\sin \beta \cdot \sinh x \sinh y} .

Using \sin \beta / \sinh z = \sin \beta_y'/ \sinh y (sine rule) once more, \sin \beta \sinh y = \sin \beta_y' \sinh z, and we obtain the clean Cho–Kim form:

\boxed{\;\; \frac{\partial \beta}{\partial x} \;=\; -\frac{\cos \beta_y'}{\sinh x}, \quad
            \frac{\partial \beta}{\partial y} \;=\; -\frac{\cos \beta_x'}{\sinh y}, \quad
            \frac{\partial \beta}{\partial z} \;=\; +\frac{1}{\sinh z} \cdot \frac{\sinh z}{\sin\beta\,\sinh x \sinh y} \cdot \sinh z \;\;}

For the third partial, \partial_z N = -\sinh z, \partial_z D = 0, so

-\sin \beta \cdot \partial_z \beta = \frac{-\sinh z}{\sinh x \sinh y}
\quad \Longrightarrow \quad
\boxed{\;\; \frac{\partial \beta}{\partial z}
   \;=\; \frac{\sinh z}{\sin \beta \cdot \sinh x \sinh y}
   \;=\; \frac{1}{\sin \beta} \cdot \frac{\sinh z}{\sinh x \sinh y} . \;\;}

Equivalently, \partial \beta / \partial z = \sinh z / (\sin \beta \sinh x \sinh y), which by the sine rule is also 1 / (\sinh x \sin \beta_y'), recovering a form symmetric to the first two.

Summary (operational form used in code). For \beta = \zeta(x, y, z) (opposite to z):

\begin{aligned}
\partial_x \beta &= \frac{1}{\sin \beta} \cdot \frac{\cosh y - \cosh x \cosh z}{\sinh^2 x \sinh y}, \\
\partial_y \beta &= \frac{1}{\sin \beta} \cdot \frac{\cosh x - \cosh y \cosh z}{\sinh x \sinh^2 y}, \\
\partial_z \beta &= \frac{1}{\sin \beta} \cdot \frac{\sinh z}{\sinh x \sinh y} .
\end{aligned}
\tag{Z}

These three lines are what the analytic kernel will evaluate (they are finite as long as no edge degenerates and \sin \beta \neq 0, which is guaranteed inside the triangle-inequality regime gated by face_angles_from_local_dofs).

---

3. Derivatives of \ell_{ij} (edge-length building blocks)

The length dispatcher lij(b_i, b_j, a_{ij}, v_i, v_j) in hyper_ideal_geometry.hpp selects among \zeta_{13}, \zeta_{14}, \zeta_{15} based on which vertices are hyper-ideal. We treat each branch separately.

3.1 Both vertices hyper-ideal: \ell_{ij} = \zeta_{13}(b_i, b_j, a_{ij})

\zeta_{13}(x, y, z) = \operatorname{arcosh}\!\left( \frac{\cosh x \cosh y + \cosh z}{\sinh x \sinh y} \right) .

Let L = \zeta_{13}, \Phi := \cosh L = (\cosh x \cosh y + \cosh z)/(\sinh x \sinh y). Then \sinh L \cdot \partial L = \partial \Phi. With

\partial_x \Phi = \frac{\sinh x \cosh y \cdot \sinh x \sinh y - (\cosh x \cosh y + \cosh z) \cosh x \sinh y}{\sinh^2 x \sinh^2 y}
= \frac{-\cosh x \cosh z - \cosh y}{\sinh^2 x \sinh y},

(using \sinh^2 x - \cosh^2 x = -1), and analogously for y, z:

\boxed{\;\;
\begin{aligned}
\partial_x \zeta_{13} &= -\frac{1}{\sinh L} \cdot \frac{\cosh y + \cosh x \cosh z}{\sinh^2 x \sinh y}, \\
\partial_y \zeta_{13} &= -\frac{1}{\sinh L} \cdot \frac{\cosh x + \cosh y \cosh z}{\sinh x \sinh^2 y}, \\
\partial_z \zeta_{13} &= +\frac{1}{\sinh L} \cdot \frac{\sinh z}{\sinh x \sinh y} .
\end{aligned}
\;\;}
\tag{Z13}

Compare (Z) and (Z13): they differ only in (i) the sign in front of \cosh z in the numerator of the first two partials, and (ii) the prefactor (1/\sin \beta vs 1/\sinh L) — which mirrors the \arccos / \operatorname{arcosh} duality. This is the dual hyperbolic law of cosines for a right-angled hexagon (Buser 1992, Geometry and Spectra of Compact Riemann Surfaces, Thm 2.4.1).

3.2 One ideal vertex: \ell_{ij} = \zeta_{14}(a_{ij}, b_j)

If vertex i is ideal (i.e. only v_j is hyper-ideal),

\zeta_{14}(x, y) = \operatorname{arcosh}\!\left( \frac{e^x + \cosh y}{\sinh y} \right) .

Let L = \zeta_{14}, \Phi = \cosh L = (e^x + \cosh y)/\sinh y. Then

\partial_x \Phi = \frac{e^x}{\sinh y},
\qquad
\partial_y \Phi = \frac{\sinh y \cdot \sinh y - (e^x + \cosh y)\cosh y}{\sinh^2 y}
       = \frac{-1 - e^x \cosh y}{\sinh^2 y} .

Hence

\boxed{\;\;
\partial_x \zeta_{14} \;=\; \frac{e^x}{\sinh L \cdot \sinh y},
\qquad
\partial_y \zeta_{14} \;=\; -\frac{1 + e^x \cosh y}{\sinh L \cdot \sinh^2 y} .
\;\;}
\tag{Z14}

Note that lij invokes \zeta_{14} with argument order (edge, vertex), so when v_i is the ideal one the call is zeta14(a_{ij}, b_j) and we must remember that \partial / \partial a_{ij} = \partial_x \zeta_{14}, \partial / \partial b_j = \partial_y \zeta_{14}. The symmetric case (when v_j is the ideal one) flips the roles of b_i / b_j.

3.3 Both vertices ideal: \ell_{ij} = \zeta_{15}(a_{ij})

\zeta_{15}(x) = 2\, \operatorname{arsinh}\!\big(e^{x/2}\big), \qquad
\partial_x \zeta_{15} = \frac{2 \cdot \tfrac{1}{2} e^{x/2}}{\sqrt{1 + e^x}}
                       = \frac{e^{x/2}}{\sqrt{1 + e^x}} .

A more numerically symmetric form uses \cosh(\zeta_{15}/2) = \sqrt{1 + e^x}, \sinh(\zeta_{15}/2) = e^{x/2}, so

\boxed{\;\;
\partial_x \zeta_{15} \;=\; \tanh\!\big(\zeta_{15}(x)/2\big) . \;\;}
\tag{Z15}

This is the cleanest form for cross-checks against the FD reference.

3.4 Combined edge-length partials

We package the per-edge \partial \ell_{ij} as a length differential

d\ell_{ij} \;=\; L^b_{ij,i}\, db_i \;+\; L^b_{ij,j}\, db_j \;+\; L^a_{ij}\, da_{ij}

with case-by-case coefficients:

Case (v_i, v_j)	L^b_{ij,i}	L^b_{ij,j}	L^a_{ij}
(hyp, hyp)	\partial_x \zeta_{13} at (b_i, b_j, a_{ij})	\partial_y \zeta_{13}	\partial_z \zeta_{13}
(ideal, hyp)	—	\partial_y \zeta_{14} at (a_{ij}, b_j)	\partial_x \zeta_{14}
(hyp, ideal)	\partial_y \zeta_{14} at (a_{ij}, b_i)	—	\partial_x \zeta_{14}
(ideal, ideal)	—	—	\partial_x \zeta_{15} at a_{ij}

A "—" entry means the corresponding DOF does not exist; the partial is identically zero on the constraint surface.

---

4. Chain-rule assembly (interior angles and dihedrals)

4.1 Interior angles \beta_i

On face f = (1, 2, 3) with edges \ell_{12}, \ell_{23}, \ell_{31}, the code computes

\beta_1 = \zeta(\ell_{12}, \ell_{31}, \ell_{23}),\quad
\beta_2 = \zeta(\ell_{23}, \ell_{12}, \ell_{31}),\quad
\beta_3 = \zeta(\ell_{31}, \ell_{23}, \ell_{12}) .

With (Z) we get, for any DOF \xi \in \{b_1, b_2, b_3, a_{12}, a_{23}, a_{31}\},

\frac{\partial \beta_1}{\partial \xi}
   \;=\; \zeta_x(\ell_{12}, \ell_{31}, \ell_{23})\, \partial_\xi \ell_{12}
   \;+\; \zeta_y(\ell_{12}, \ell_{31}, \ell_{23})\, \partial_\xi \ell_{31}
   \;+\; \zeta_z(\ell_{12}, \ell_{31}, \ell_{23})\, \partial_\xi \ell_{23} ,
\tag{B1}

and cyclically for \beta_2, \beta_3. Each \partial_\xi \ell_{e} is read from the table in §3.4: of the six DOFs only the three that touch edge e contribute (i.e. b_{e^-}, b_{e^+}, a_e).

4.2 Dihedral angles \alpha_{ij} — general structure

alpha_ij(...) has four branches depending on (v_i, v_j, v_k). Define the "$\sigma$-triangle" attached to vertex i:

s_i \;=\; \sigma_i(a_{ij}, a_{ki}, a_{jk}; v_j, v_k), \qquad
s_{ij} \;=\; \sigma_{ij}(a_{ij}, b_i, b_j; v_j), \qquad
s_{ik} \;=\; \sigma_{ij}(a_{ki}, b_i, b_k; v_k) .

Then in the hyper-ideal $v_i$ branch

\alpha_{ij} \;=\; \zeta(s_i, s_{ij}, s_{ik})
\tag{A.v_i}

— note that the dihedral \alpha_{ij} is computed as an interior angle in the half-triangle at vertex i, with s_i playing the role of side x, s_{ij} of y and s_{ik} of z.

The branches v_j (hyper-ideal but v_i ideal), v_k (one level of recursion), and all ideal (closed form) follow the same pattern.

4.3 Derivatives of \sigma_i and \sigma_{ij}

Both are themselves dispatchers over \zeta_{13}, \zeta_{14}, \zeta_{15}.

\sigma_{ij}(a_{ij}, b_i, b_j; v_j).

• If v_j hyper-ideal: \sigma_{ij} = \zeta_{13}(a_{ij}, b_i, b_j) with partials (\partial_x \zeta_{13}, \partial_y \zeta_{13}, \partial_z \zeta_{13}) evaluated at (a_{ij}, b_i, b_j). Note the argument order: the first slot of \zeta_{13} is a_{ij}, not b_i.
• If v_j ideal: \sigma_{ij} = \zeta_{14}(-a_{ij}, b_i). Then \partial_{a_{ij}} \sigma_{ij} = -\partial_x \zeta_{14}(-a_{ij}, b_i), \partial_{b_i} \sigma_{ij} = \partial_y \zeta_{14}(-a_{ij}, b_i), \partial_{b_j} \sigma_{ij} = 0.

\sigma_i(a_{ij}, a_{ki}, a_{jk}; v_j, v_k).

• (v_j, v_k) = (hyp, hyp): \sigma_i = \zeta_{13}(a_{ij}, a_{ki}, a_{jk}). Partials direct.
• (v_j, v_k) = (hyp, ideal): \sigma_i = \zeta_{14}(a_{jk} - a_{ki}, a_{ij}). Then \partial_{a_{jk}} = \partial_x \zeta_{14}, \partial_{a_{ki}} = -\partial_x \zeta_{14}, \partial_{a_{ij}} = \partial_y \zeta_{14}.
• (v_j, v_k) = (ideal, hyp): \sigma_i = \zeta_{14}(a_{jk} - a_{ij}, a_{ki}). Sign pattern is the mirror image.
• (v_j, v_k) = (ideal, ideal): \sigma_i = \zeta_{15}(a_{jk} - a_{ij} - a_{ki}). Partials are (\partial_x \zeta_{15}, -\partial_x \zeta_{15}, -\partial_x \zeta_{15}) along (a_{jk}, a_{ij}, a_{ki}).

These nine sub-branches are the bulk of the per-face symbolic work.

4.4 Chain rule for \alpha_{ij} in the $v_i$-branch

By (A.v_i) and (Z) evaluated at (s_i, s_{ij}, s_{ik}),

d\alpha_{ij} \;=\;
\zeta_x|_{(s_i, s_{ij}, s_{ik})} ds_i
\;+\; \zeta_y|_{(s_i, s_{ij}, s_{ik})} ds_{ij}
\;+\; \zeta_z|_{(s_i, s_{ij}, s_{ik})} ds_{ik} ,
\tag{A1}

and each ds_i, ds_{ij}, ds_{ik} expands via §4.3. The result is a linear combination of the six face DOFs with closed-form coefficients.

For the $v_j$-branch swap the roles (i \leftrightarrow j); the formulas are identical up to a permutation of $\sigma$-indices.

4.5 The recursive $v_k$-branch

When v_i, v_j are both ideal but v_k is hyper-ideal, the code returns

\alpha_{ij} \;=\; \pi - \alpha_{jk}' - \beta_j ,

where \alpha_{jk}' is a recursive call into the $v_j$-branch (since the cyclic role rotation \, (i, j, k) \mapsto (j, k, i) makes the second- position vertex v_k, which is hyper-ideal, trigger the $v_i$-branch on the recursion). Differentiating,

d\alpha_{ij} \;=\; -\, d\alpha_{jk}' \;-\; d\beta_j ,
\tag{A2}

so we obtain d\alpha_{ij} by computing d\alpha_{jk}' via (A1) (with the rotated DOF identification) and subtracting d\beta_j from §4.1. This single-level recursion is finite because the recursive call descends into the $v_i$-branch (whose v_i is now the original v_k, which is hyper-ideal by hypothesis); see the comment "never more than one level deep" in hyper_ideal_geometry.hpp line 106.

4.6 The all-ideal branch (closed form)

When all three vertices are ideal, \alpha_{ij} = \tfrac{1}{2}(\pi + \beta_k - \beta_i - \beta_j). Hence

d\alpha_{ij} \;=\; \tfrac{1}{2}(d\beta_k - d\beta_i - d\beta_j) ,
\tag{A3}

with each d\beta_\bullet from (B1). In this case the DOF vector collapses to (a_{12}, a_{23}, a_{31}) (the $b$'s do not exist), and the six-by-six block reduces to a non-trivial 3 \times 3 block embedded along the $a$-axes.

---

5. The per-face 6 \times 6 block

5.1 Layout

Order the inputs and outputs of face_angles_from_local_dofs as

x_{\text{loc}} \;=\; (b_1, b_2, b_3, a_{12}, a_{23}, a_{31})^\top,
\qquad
y_{\text{loc}} \;=\; (\beta_1, \beta_2, \beta_3, \alpha_{12}, \alpha_{23}, \alpha_{31})^\top .

The face Jacobian is

J(f) \;=\; \frac{\partial y_{\text{loc}}}{\partial x_{\text{loc}}}
        \;=\;
\begin{pmatrix}
J^{\beta b} & J^{\beta a} \\
J^{\alpha b} & J^{\alpha a}
\end{pmatrix} \in \mathbb{R}^{6 \times 6} .

The 3 \times 3 sub-blocks are:

• J^{\beta b}_{ij} = \partial \beta_i / \partial b_j, computed from (B1) using the L^b coefficients from §3.4. Sparse: only the two edges (i, \cdot) adjacent to vertex i contribute the b_j derivative (via \partial \ell_{ij} / \partial b_j), so J^{\beta b}_{ij} has at most two nonzero edge-length terms per (i, j) pair.
• J^{\beta a}_{ie} = \partial \beta_i / \partial a_e, also from (B1): exactly one $\zeta$-derivative slot per (i, e) pair, since each a_e affects exactly one edge length and each \ell_e enters \beta_i in exactly one slot.
• J^{\alpha b}_{e j} = \partial \alpha_e / \partial b_j: computed from the branch-appropriate identity among (A1), (A2), (A3). For (A1), only s_{ij} and s_{ik} depend on $b$'s, so the row collapses to \zeta_y \cdot \partial_{b_j} s_{ij} + \zeta_z \cdot \partial_{b_j} s_{ik}.
• J^{\alpha a}_{e e'} = \partial \alpha_e / \partial a_{e'}: this is the densest block; both s_i and s_{ij}, s_{ik} depend on $a$'s, so all three $\zeta$-slots contribute.

5.2 Closed-form formulas, hyper-ideal triangle (all v_i variable)

Let \beta_1 = \zeta(\ell_{12}, \ell_{31}, \ell_{23}) (and cyclic), and write the $\zeta$-partials at vertex 1 as \zeta_x^{(1)} := \zeta_x(\ell_{12}, \ell_{31}, \ell_{23}), etc.; and the length partials from (Z13) as L^b_{e, i}, L^a_e etc. Then explicitly:

\begin{aligned}
\frac{\partial \beta_1}{\partial b_1} &=
   \zeta_x^{(1)} L^b_{12,1} \;+\; \zeta_y^{(1)} L^b_{31,1} ,
   &\quad
\frac{\partial \beta_1}{\partial b_2} &= \zeta_x^{(1)} L^b_{12,2} , \\
\frac{\partial \beta_1}{\partial b_3} &= \zeta_y^{(1)} L^b_{31,3} ,
   &\quad
\frac{\partial \beta_1}{\partial a_{12}} &= \zeta_x^{(1)} L^a_{12} , \\
\frac{\partial \beta_1}{\partial a_{31}} &= \zeta_y^{(1)} L^a_{31} ,
   &\quad
\frac{\partial \beta_1}{\partial a_{23}} &= \zeta_z^{(1)} L^a_{23} .
\end{aligned}
\tag{Bblock}

(Note: \beta_1 does not depend on a_{23} via the $b$'s — only the direct $\ell_{23}$-dependence through \zeta_z^{(1)} — so the last line is the only $a_{23}$-coupling of \beta_1.)

The full J^{\beta b} row for vertex 1 is therefore:

J^{\beta b}_{1,*} \;=\;
\big(\;
   \zeta_x^{(1)} L^b_{12,1} + \zeta_y^{(1)} L^b_{31,1} ,\;\;
   \zeta_x^{(1)} L^b_{12,2} ,\;\;
   \zeta_y^{(1)} L^b_{31,3}
\;\big) ,

and similarly for rows 2 and 3 by cyclic permutation.

5.3 Closed-form formulas, dihedral row in $v_1$-branch

Let \beta := \alpha_{12} = \zeta(s_1, s_{12}, s_{13}) (in the v_1 branch, the i of \sigma is vertex 1, and the two edges through vertex 1 are 12 and 31 = 13). With \zeta_x^{(\alpha_{12})} etc. denoting the three partials of \zeta evaluated at (s_1, s_{12}, s_{13}):

\frac{\partial \alpha_{12}}{\partial \xi}
   \;=\; \zeta_x^{(\alpha_{12})} \frac{\partial s_1}{\partial \xi}
   \;+\; \zeta_y^{(\alpha_{12})} \frac{\partial s_{12}}{\partial \xi}
   \;+\; \zeta_z^{(\alpha_{12})} \frac{\partial s_{13}}{\partial \xi} ,

where \xi ranges over the six face DOFs. Substituting s_1 = \zeta_{13}(a_{12}, a_{31}, a_{23}) (assuming all v_j, v_k hyper-ideal so we are in the (\text{hyp, hyp}) sub-branch of \sigma_1), s_{12} = \zeta_{13}(a_{12}, b_1, b_2), s_{13} = \zeta_{13}(a_{31}, b_1, b_3), the six entries are:

\begin{array}{l|l}
\xi & \partial \alpha_{12} / \partial \xi \\\hline
b_1 & \zeta_y^{(\alpha_{12})} \partial_y \zeta_{13}|_{s_{12}} + \zeta_z^{(\alpha_{12})} \partial_y \zeta_{13}|_{s_{13}} \\
b_2 & \zeta_y^{(\alpha_{12})} \partial_z \zeta_{13}|_{s_{12}} \\
b_3 & \zeta_z^{(\alpha_{12})} \partial_z \zeta_{13}|_{s_{13}} \\
a_{12} & \zeta_x^{(\alpha_{12})} \partial_x \zeta_{13}|_{s_1} + \zeta_y^{(\alpha_{12})} \partial_x \zeta_{13}|_{s_{12}} \\
a_{31} & \zeta_x^{(\alpha_{12})} \partial_y \zeta_{13}|_{s_1} + \zeta_z^{(\alpha_{12})} \partial_x \zeta_{13}|_{s_{13}} \\
a_{23} & \zeta_x^{(\alpha_{12})} \partial_z \zeta_{13}|_{s_1}
\end{array}
\tag{Ablock}

The rows for \alpha_{23} and \alpha_{31} follow by cyclic permutation (1, 2, 3) \mapsto (2, 3, 1) \mapsto (3, 1, 2).

5.4 Symmetry check (Schläfli)

By (S2), the 6 \times 6 block H(f) = J(f) — interpreted as a Hessian-of-energy block via the Springborn identification \beta_i = \partial U(f)/\partial b_i, \alpha_e = \partial U(f)/\partial a_e — satisfies J(f) = J(f)^\top. Equivalently the four sub-blocks obey

J^{\beta b} = (J^{\beta b})^\top, \quad J^{\alpha a} = (J^{\alpha a})^\top,
\quad J^{\alpha b} = (J^{\beta a})^\top .
\tag{Sym}

This is a powerful numerical sanity check: any analytic formula that violates (Sym) by more than rounding error has a derivation error somewhere in §3–§4.

The off-diagonal cross-check \partial \alpha_{12} / \partial b_3 \;\stackrel{!}{=}\; \partial \beta_3 / \partial a_{12} is particularly instructive: the left-hand side is computed via the $\sigma$-triangle at vertex 1 (last row of Ablock-style table for \alpha_{12}), the right-hand side via the auxiliary triangle at vertex 3. The two routes coincide only because of (S2).

5.5 Scatter to global Hessian

Once J(f) is built, the scatter step is identical to the existing hyper_ideal_hessian_block_fd: for each (i, j) \in \{1, \dots, 6\}^2, add J_{ij}(f) to H[\text{glb}(i), \text{glb}(j)], where \text{glb} maps the local face slot to its global DOF index via v_idx / e_idx. Pinned slots (-1) are skipped, exactly as in the existing block-FD implementation.

---

6. Acceptance criteria

These mirror doc/roadmap/research-track.md Phase 9b-analytic:

1. Per-case derivative cross-checks against block-FD. For each of the four building blocks $(\partial \beta / \partial \ell, \partial \ell / \partial \xi, \partial \sigma / \partial \xi, \partial \alpha / \partial \xi)$, sample 100 random DOF vectors in

   • all-hyper-ideal regime (v_1, v_2, v_3 variable),
   • one-ideal regime (each of three rotations of v_i pinned),
   • two-ideal regime (each of three rotations of two v_i pinned). Require \max_{ij} |J_{ij}^{\text{analytic}} - J_{ij}^{\text{FD}}| \le 10^{-6} with central-difference step \varepsilon = 10^{-5}.

2. Schläfli symmetry (gauge-free). For each face f at random x, verify \|J(f) - J(f)^\top\|_\infty \le 10^{-10}. Symbolically this is automatic by (S2); numerically it pins down sign/branch bugs.

3. Gauge null space. Let \mathbf{1}_b be the constant-b Möbius dilation mode (i.e. the vector with 1 in every $b$-slot and 0 in every $a$-slot). Require \|H \mathbf{1}_b\|_\infty \le 10^{-10} on the all-hyper-ideal mesh after assembly (this is the only Möbius mode that survives in the closed-surface case).

4. PSD on the interior. For random DOFs inside the convex domain (Springborn 2020 §4.3), verify \lambda_{\min}(H) \ge -10^{-12} (numerical zero). Closed-form: E is convex on the admissible cone by Springborn 2020 Theorem 4.2.

5. Speed-up. Wall-clock benchmark on tetrahedron, octahedron, icosahedron, and a 200-vertex genus-2 test mesh: analytic Hessian assembly \ge 3\times faster than hyper_ideal_hessian_block_fd, with target \sim 6\times asymptotically (one face costs O(6 \cdot 36) FD evaluations of face_angles_from_local_dofs currently; analytic costs \sim 100 scalar transcendentals per face).

---

7. Implementation outline

Header: code/include/hyper_ideal_hessian_analytic.hpp (new).

// hyper_ideal_hessian_analytic.hpp  — Phase 9b-analytic
//
// Closed-form 6x6 per-face Hessian via Schläfli + chain rule.
// See doc/math/hyperideal-hessian-derivation.md for the derivation.

#include "hyper_ideal_geometry.hpp"
#include "hyper_ideal_functional.hpp"

namespace conformallab {

// Building-block derivatives (Section 3 of the note).
struct ZetaPartials   { double dx, dy, dz; };  // d zeta(x, y, z)
struct Zeta13Partials { double dx, dy, dz; };  // d zeta13
struct Zeta14Partials { double dx, dy;     };  // d zeta14(x, y)

inline ZetaPartials   d_zeta  (double x, double y, double z, double beta);
inline Zeta13Partials d_zeta13(double x, double y, double z, double L);
inline Zeta14Partials d_zeta14(double x, double y, double L);
inline double         d_zeta15(double x);             // returns tanh(L/2)

// Edge-length differential coefficients (table in 3.4).
struct EdgeLenDiff {
    double dbi;    // d l_ij / d b_i   (zero if v_i ideal)
    double dbj;    // d l_ij / d b_j   (zero if v_j ideal)
    double daij;   // d l_ij / d a_ij
};
EdgeLenDiff d_lij(double bi, double bj, double aij, bool vi, bool vj);

// 6x6 face Jacobian: rows (beta1, beta2, beta3, alpha12, alpha23, alpha31),
// cols (b1, b2, b3, a12, a23, a31).  Same input contract as
// face_angles_from_local_dofs.
struct FaceJacobian6 { double M[6][6]; };

FaceJacobian6 analytic_face_jacobian(
    double b1, double b2, double b3,
    double a12, double a23, double a31,
    bool v1b, bool v2b, bool v3b);

// Drop-in replacement for hyper_ideal_hessian_block_fd.
// Same scatter loop; only the inner block changes from FD to closed form.
Eigen::SparseMatrix<double> hyper_ideal_hessian_analytic(
    const ConformalMesh& mesh,
    const std::vector<double>& x,
    const HyperIdealMaps& m);

} // namespace conformallab

The body of analytic_face_jacobian follows the case dispatch:

1. clamp inputs same as face_angles_from_local_dofs
2. compute l12, l23, l31 via lij; bail out on triangle-inequality break
3. compute beta1, beta2, beta3 via zeta
4. for each edge e, compute EdgeLenDiff via d_lij(...)
5. for each vertex i, fill row i of J^{beta b} and J^{beta a}
   using (B1) with d_zeta(...) evaluated at the angle of vertex i
6. for each edge e, dispatch on (v_i, v_j, v_k) and compute alpha_e row:
   - (v_i hyp branch): use (A.v_i), expand s_i, s_ij, s_ik via 4.3
   - (v_j hyp branch): symmetric
   - (v_k recursive branch): compute alpha_jk row first, then
     subtract beta_j row and negate (A2)
   - (all-ideal branch): combine three beta rows via (A3)
7. assert |J - J^T|_inf < 1e-10  (DEBUG only)
8. return J

The global assembler then reuses the existing scatter loop from hyper_ideal_hessian_block_fd; only the inner per-face block construction changes.

---

8. References

• Schläfli, L. (1858/60). On the multiple integral \int dx\, dy \cdots dz whose limits are p_1 = a_1 x + b_1 y + \cdots + h_1 z > 0, p_2 > 0, \dots, p_n > 0 and x^2 + y^2 + \cdots + z^2 < 1. Quarterly Journal of Pure and Applied Mathematics 2, 269–301; 3, 54–68, 97–108. Reprinted in Gesammelte Mathematische Abhandlungen, Band I, Birkhäuser 1950. → first-order Schläfli identity, equation (S1) above.

• Milnor, J. (1982). Hyperbolic geometry: the first 150 years. Bulletin AMS 6, 9–24. → modern restatement of (S1) and discussion of its role in hyperbolic volume.

• Vinberg, E. B. (ed.) (1993). Geometry II: Spaces of constant curvature. Encyclopaedia of Mathematical Sciences vol. 29, Springer. Ch. 7 covers the polyhedral Schläfli identity and its derivation via the Gauss–Bonnet formula for hyperbolic polytopes.

• Cho, Y. & Kim, H. (1999). On the volume formula for hyperbolic tetrahedra. Discrete & Computational Geometry 22, 347–366. → explicit \partial \alpha / \partial a, \partial \beta / \partial b, etc. for hyperbolic tetrahedra (used in §2.2 above).

• Rivin, I. (1994). Euclidean structures on simplicial surfaces and hyperbolic volume. Annals of Math. 139, 553–580. → second-order Schläfli (S2) in the convex-polyhedron setting.

• Buser, P. (1992). Geometry and Spectra of Compact Riemann Surfaces. Birkhäuser. Thm 2.4.1 — right-angled hexagon identities (used in §3.1 to interpret \zeta_{13}).

• Glickenstein, D. (2011). Discrete conformal variations and scalar curvature on piecewise flat manifolds. J. Diff. Geom. 87, 201–238. Equation (4.6) and §4 — cone-vertex variant of the angle-length derivative formulas (mirror Hessian for inversive distance).

• Springborn, B. (2020). Ideal Hyperbolic Polyhedra and Discrete Uniformization. Discrete & Computational Geometry 64, 63–108. §4 — hyper-ideal energy and its gradient/Hessian; §4.3 — convexity (PSD criterion).

• Bobenko, A. I., Pinkall, U. & Springborn, B. (2015). Discrete conformal maps and ideal hyperbolic polyhedra. Geometry & Topology 19, 2155–2215. — earlier exposition of the truncated-tetrahedron construction that face_angles_from_local_dofs realises.

---

Appendix A. Sign conventions and pitfalls

A few places where the derivation interacts subtly with the code:

• zeta(l_jk, l_ki, l_ij) returns the interior angle at the vertex between sides l_{jk} and $l_{ki}$, hence opposite to l_{ij}. So when reading off \zeta_x, \zeta_y, \zeta_z for $\beta_1 = \zeta(\ell_{12}, \ell_{31}, \ell_{23})$, the slot x matches \ell_{12}, y matches \ell_{31}, z matches \ell_{23}. The Cho–Kim formulas (Z) must be applied with exactly this slot mapping.

• sigma_ij(a_{ij}, b_i, b_j, v_j) — when v_j is ideal, the call is zeta14(-a_{ij}, b_i), so the partial w.r.t. a_{ij} carries a negative sign. This is a frequent source of off-by-sign bugs.

• The recursive $v_k$-branch in alpha_ij returns \pi - \alpha_{jk}' - \beta_j, not \pi - \alpha_{jk}' - \beta_k. This reflects the geometry of the truncated tetrahedron and is not the cyclically-symmetric formula one might guess; see Springborn 2020 Fig. 3.

• The "all-ideal" closed form \alpha_{ij} = \tfrac{1}{2}(\pi + \beta_k - \beta_i - \beta_j) is dual: \alpha + \beta on each ideal vertex always sums to \pi/2 for a planar triangle, so the dihedral degenerates to the half-angle identity. Differentiating loses any explicit \ell dependence, hence (A3) is the simplest of the four branches.

Appendix B. Numerical guard zones

The block-FD currently used as reference clamps:

• b_i < 0 \;\to\; 0.01 on variable vertices,
• a_e < 0 \;\to\; 0 on edges between two variable vertices.

The analytic kernel must reproduce these clamps before computing any derivatives; otherwise the FD-vs-analytic cross-check (Acceptance §1) will fail at boundary inputs. Inside the clamped region, all \sinh \ell, \sin \beta denominators in (Z), (Z13), (Z14), (Z15) are bounded away from zero, so the closed-form expressions are numerically stable. The degenerate triangle-inequality branches in face_angles_from_local_dofs (lines 174–183) set angles to 0 or \pi, which makes some $\zeta$-partials infinite; in those branches the Jacobian is undefined and we return the all-zero block (matching the non-smooth boundary of the domain). This is also the policy of the block-FD reference, so the cross-check stays consistent.