# Puiseux Arclength Continuation for Implicit Layers

- ID: 473
- Canonical URL: https://synthcore.org/idea/473/puiseux-arclength-continuation-for-implicit-layers
- API JSON: https://synthcore.org/api/idea/473.json
- API Markdown: https://synthcore.org/api/idea/473.md
- Verification status: mechanism_failed
- Source: [arXiv:2607.06329](https://arxiv.org/abs/2607.06329)
- Category: dynamics
- Solves: stability, speedup, scalability
- ML areas: training-dynamics, optimizer, ssm
- Math tags: dynamical-systems, numerical-analysis, approximation-theory, linear-algebra
- Ratings: usefulness 6/10; difficulty 6/10; novelty 8/10

## Idea description

Replace the usual linear predictor in continuation of an implicit neural state with a fractional-power predictor fitted from recent states, then correct the prediction using a pseudo-arclength constraint. This is designed for equilibrium layers, implicit sequence models, or homotopy training schedules where the state Jacobian becomes nearly singular and ordinary Newton correction or fixed-point iteration becomes unstable.

## Mathematical statement

The paper starts from a homotopy curve defined by h(z(t),t)=0, where z is the implicit state and t is a continuation parameter, and differentiates it to obtain the Davidenko ODE dot z(t)=-h_z(z(t),t)^(-1)h_t(z(t),t). This formula explains why direct tangent integration fails when h_z is singular. Near a singular endpoint, the branch is represented by a Puiseux expansion z(t)=a_0+a_1 t^(k_1/c)+a_2 t^(k_2/c)+a_3 t^(k_3/c)+..., where a_i are state-space coefficient vectors, c is a positive integer, and k_i are increasing positive integers. The exponents need not be 1/c,2/c,3/c because they belong to a value semigroup. Equation (22) estimates coefficients from state and derivative observations: the matrix with rows [1, lambda^(k_1/c), lambda^(k_2/c), lambda^(k_3/c)] and [0, (k_1/c)lambda^(k_1/c), (k_2/c)lambda^(k_2/c), (k_3/c)lambda^(k_3/c)] multiplies [a_0, a_1 t_1^(k_1/c), a_2 t_1^(k_2/c), a_3 t_1^(k_3/c)]^T and equals [z_1, t_1 dot z_1, z_2, t_2 dot z_2]^T, with lambda=t_2/t_1. The neural adaptation uses fitted exponents alpha_i=k_i/c and predicts z_pred(t)=a_0+sum_i a_i t^(alpha_i), while replacing exact tangents with secants or Jacobian-vector solves. A pseudo-arclength correction augments the residual equation with a tangent constraint, allowing continuation even when h_z is rank-deficient.

## Key formulas

- $$h(z(t),t)=0,\qquad \dot{z}(t)=-h_z(z(t),t)^{-1}h_t(z(t),t).$$
- $$z(t)=a_0+a_1t^{k_1/c}+a_2t^{k_2/c}+a_3t^{k_3/c}+\cdots,$$
- $$\begin{bmatrix}1&\lambda^{k_1/c}&\lambda^{k_2/c}&\lambda^{k_3/c}\\0&\frac{k_1}{c}\lambda^{k_1/c}&\frac{k_2}{c}\lambda^{k_2/c}&\frac{k_3}{c}\lambda^{k_3/c}\end{bmatrix}\begin{bmatrix}a_0\\a_1t_1^{k_1/c}\\a_2t_1^{k_2/c}\\a_3t_1^{k_3/c}\end{bmatrix}=\begin{bmatrix}z_1\\t_1\dot z_1\\z_2\\t_2\dot z_2\end{bmatrix},\qquad \lambda=\frac{t_2}{t_1}.$
- $$\begin{bmatrix}h_z(z,t)&h_t(z,t)\\q_z^{\mathsf T}&q_t\end{bmatrix}\begin{bmatrix}\Delta z\\\Delta t\end{bmatrix}=-\begin{bmatrix}h(z,t)\\q_z^{\mathsf T}(z-z_{\mathrm{pred}})+q_t(t-t_{\mathrm{pred}})\end{bmatrix},$$ where $(q_z,q_t)$ is the normalized secant tangent from the previous two continuation points.

## Implementation notes

Integrate this at the implicit-state solve of an equilibrium layer, for example z=phi_theta(z,x,t) with residual h(z,t)=z-phi_theta(z,x,t); t interpolates a stable initialization or regularized problem at t=1 toward the target problem at small t. Maintain the last four accepted pairs (t_j,z_j). Before each step, estimate a local exponent by fitting log(||z_j-z_{j-1}||+epsilon)=b+alpha log(t_j+epsilon) with least squares over the recent differences. Clamp alpha to [0.1,4], reject unstable fits, and retain up to three distinct exponents alpha_i. Solve the displayed coefficient system using recent states and finite-difference tangents dot z_j=(z_j-z_{j-1})/(t_j-t_{j-1}); with more history than equations, use regularized least squares. Predict at t_new=t_j-delta_t using z_pred=a_0+sum_i a_i t_new^(alpha_i). Correct jointly in (z,t) with Newton or damped least squares on h(z,t)=0 plus the pseudo-arclength equation q_z^T(z-z_pred)+q_t(t-t_new)=0. Compute h_z and h_t with automatic differentiation, but estimate exponents and coefficients empirically; the paper contributes the fractional-power model and derivative scaling rather than neural-specific hyperparameters. Use step rejection when Newton fails, and adapt delta_t by doubling after success or halving after failure. First test a small 2- or 3-layer DEQ on MNIST or CIFAR-10, comparing this method with linear prediction plus fixed-parameter Newton over a homotopy from t=1 to 10^-4. Record failed solves, Newton iterations, accepted step sizes, wall-clock time, final residual, and task accuracy. A positive result is at least 1.5x larger average continuation steps near rank deficiency and fewer failed corrections at equal residual tolerance.

## Verification

- Status: mechanism_failed
- Mechanism evidence: yes
- Mechanism confirmed: no
- Verdict: Built a scalar Puiseux/arclength continuation MVP for the singular implicit branch h(z,t)=z^2-t, including coefficient fitting, exponent estimation, pseudo-arclength Newton correction, and a linear-predictor baseline. The coefficient-system check passed with 6.7e-16 equation error and 7.2e-16 scaled-coefficient error. On 80 matched geometric steps from t=1 to 1e-8, Puiseux correction averaged 2.71 Newton iterations versus 3.26 for linear prediction, but both had zero failures; the idea was ~31x slower in this tiny implementation, prediction error was slightly worse, and exponent fitting degraded near machine precision, so there is no clear robust promised effect.

### Mechanism check

- Verdict: Built a scalar Puiseux/arclength continuation MVP for the singular implicit branch h(z,t)=z^2-t, including coefficient fitting, exponent estimation, pseudo-arclength Newton correction, and a linear-predictor baseline. The coefficient-system check passed with 6.7e-16 equation error and 7.2e-16 scaled-coefficient error. On 80 matched geometric steps from t=1 to 1e-8, Puiseux correction averaged 2.71 Newton iterations versus 3.26 for linear prediction, but both had zero failures; the idea was ~31x slower in this tiny implementation, prediction error was slightly worse, and exponent fitting degraded near machine precision, so there is no clear robust promised effect.
- Confidence: 8/10
- Limitations: Only a scalar analytic branch was tested, not a neural equilibrium layer, multidimensional rank-deficient Jacobian, adaptive step rejection, GPU execution, training accuracy, or wall-clock performance at realistic scale. The endpoint exponent estimate is contaminated by floating-point/finite-difference effects, and the benchmark uses a fixed schedule rather than a full adaptive continuation policy.

## Artifacts

- [puiseux_experiment.py](https://synthcore.org/code/1350/puiseux_experiment.py)
- [report.md](https://synthcore.org/code/1350/report.md)
- [results.json](https://synthcore.org/code/1350/results.json)
- [Download all files as ZIP](https://synthcore.org/download/1350)

## Disclaimer

AI-generated research hypothesis, automatically tested. Not peer-reviewed.