{"id":"4ea94541-ff3e-4893-b8cb-30eaf6a94aaa","arxiv_id":"2509.01855","paper_version":1,"verdict":"REJECT","confidence":"MODERATE","novelty_score":6.0,"correctness_risk":"high","formal_verification":"none","parameter_count":1,"one_line_summary":"A matrix-free Birkhoff-Chebyshev solver with FFT-based matrix-vector products and a near-perfect preconditioner solves a nonconvex orbit-transfer problem on 1,048,577 grid points with mesh-independent iteration counts.","lead":"The authors build a matrix-free trajectory optimization solver that combines a Birkhoff-Chebyshev discretization, FFT-based matrix-vector products, and a specially designed lower-triangular preconditioner, and they report solving a nonlinear orbit-transfer problem on 1,048,577 grid points on a laptop. If the speed claims hold, optimal control could be applied at a scale previously limited to much smaller grids.","discovery_kind":"new_method","skeptic_critique":{"model":"deepseek-v4-flash","headline":"The linear solver and preconditioner are derived only for Eq. (28); Table 1 reports KKT solves for the free-time, endpoint-constrained problem (39), whose KKT system is never given—and FastPinvY's printed update is inconsistent with P^{-1}.","rationale":"The paper's FastBV derivation is internally coherent and the condition-number/matrix-free motivation is reasonable; that is real evidence for the FFT product. But the central numerical demonstration depends on a linear solver whose target system is never specified. The reader's weakest assumption is on target: Eq. (28) is a scalar fixed-endpoint Newton system, while (39) is a free-time multi-state endpoint-constrained problem. No code, data, or absolute timings are included, and the one recurrence needed to make the preconditioner work is inconsistent as printed. These are not stylistic objections: they block reproduction of the exact 37-solve counts and the O(1) scaling. I agree with rejection of the current version; a revision that derives the full KKT system, corrects FastPinvY, and provides reproducible artifacts could change that.","tokens_in":20068,"tokens_out":8404,"duration_ms":92873,"concrete_test":"One check: write out the discrete KKT system for problem (39) and re-run the Table 1 experiment using the paper's algorithms. Before scaling, for N=10 with random diagonal ∂Xf, verify P·FastPinvY(e_j)=e_j exactly as printed; if the residual is not machine precision, FastPinvY is wrong. If the derived KKT system requires different block operations or the solve count changes, the million-point claim is unsupported.","verdict_should_be":"UNCHANGED","load_bearing_attack":"All of the matrix-free machinery—FastAX, FastPinvY, FastLinSol—is derived for Eq. (28), the Newton system [I−Ba∂Xf]δX = xa b for the distilled scalar, fixed-endpoint problem. The demonstrated orbit transfer (39) has four states, one control, free final time T, and endpoint constraints. Its KKT system would contain state/control variation blocks, adjoint variables, endpoint multipliers, and dT coupling; that system is never written down. Therefore the 'Matrix-free KKT solves' in Table 1 (37 for every N) and Figure 6 cannot be checked against the paper, and the scalar-diagonal argument for ∂Xf does not automatically extend to the 4N×4N dynamics Jacobian plus control/endpoint/free-time blocks. This is load-bearing because the O(1)-iteration and million-point claims depend on solving that full KKT system, not the simplified Eq. (28). Compounding the gap, FastPinvY as printed does not invert P = I−B̃a∂Xf: forward substitution requires s ← s + wk dk ξk, but the algorithm uses s ← s + wk dk Υk and solves with minus s, so the preconditioned operator AP^{-1} in Eq. (33) is not the one implemented.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"The paper claims a first-of-a-kind fast trajectory optimization solver capable of solving a nontrivial, nonlinear, nonconvex orbit-transfer problem over more than a million CGL grid points on an ordinary laptop. The enabling components are: (i) a Birkhoff-Chebyshev discretization, (ii) an FFT-based O(N log N) algorithm (FastBV) for the Birkhoff matrix-vector product, (iii) a lower-triangular Birkhoff-like preconditioner built from Clenshaw-Curtis weights, and (iv) a matrix-free preconditioned iterative solver (FastLinSol) for a simplified Newton system [I - Ba ∂X f]δX = xa b. The paper also reports mesh-independent SQP/KKT iteration counts (Table 1) for a 4-state, 1-control free-time orbit transfer problem, and displays a solution over N+1 = 1,048,577 points. The central advertised claims are O(N log N) per matrix-vector product, O(1) iteration counts in N, and the first million-point solution of a nonlinear nonconvex trajectory optimization problem.","tokens_in":20340,"tokens_out":5132,"duration_ms":56819,"significance":"If all claims held and were fully substantiated, the contribution would be significant: an O(N log N) matrix-free Birkhoff product with a mesh-independent preconditioner would remove a major bottleneck in pseudospectral optimal control and make million-point solutions practical. The FastBV derivation appears plausible and internally consistent, and the idea of a lower-triangular Birkhoff-like preconditioner is attractive. However, the paper as written does not deliver a verifiable version of those claims for the demonstrated problem: the linear solver is derived only for a scalar, fixed-endpoint system, while the reported KKT-solve counts are for a free-time, endpoint-constrained multi-state problem whose KKT system is never written down. In addition, the printed FastPinvY algorithm does not compute the inverse of the advertised preconditioner. These are load-bearing gaps in the central claim. The numerical demonstration would be impressive if the missing derivations, corrected algorithm, and verification metrics were supplied and ideally accompanied by reproducible code.","major_comments":[{"comment":"The printed FastPinvY algorithm is inconsistent with the preconditioner P := I - Btilde^a ∂X f defined in Eq. (38). Since Btilde^a is lower triangular with diagonal entries w_k/2 and subdiagonal entries w_j, forward substitution for P ξ = Υ requires s_k = Σ_{j<k} w_j d_j ξ_j and ξ_k = (1 - w_k d_k/2)^{-1}(Υ_k + s_k), followed by s_{k+1} = s_k + w_k d_k ξ_k. The algorithm as printed uses (Υ_k - s) and updates s with w_k (∂x_k f) Υ_k, not with w_k d_k ξ_k. It therefore does not compute P^{-1}Υ. This makes Eq. (33) and all subsequently reported preconditioned iteration counts unverifiable. Please correct the algorithm and, if possible, provide a small reproducible test confirming that P^{-1} is computed correctly.","section":"FastPinvY (Section 'A Preconditioned Algorithm for Solving Eq. (28)')"},{"comment":"The KKT system for the demonstrated orbit-transfer problem (39) is never stated. Problem (39) has four states, one control, free final time T, and endpoint constraints; its KKT system will contain state/control variation blocks, adjoint variables, endpoint multipliers, and dT coupling. Table 1 reports 'Matrix-free KKT solves' with counts 37 for every N, but the system being solved is not defined, and no algorithm is given for applying the FastAX/FastPinvY/FastLinSol machinery to that system. The scalar diagonal argument for ∂X f in Eq. (27) does not automatically extend to the 4N×4N dynamics Jacobian plus control/endpoint/free-time blocks. The O(1) iteration claim and the million-point result therefore rest on an unspecified computation. The authors should write out the KKT system, specify the SQP iteration, and report per-iteration linear-solve counts for that system.","section":"Section 'An Illustrative Nonlinear Nonconvex Astrodynamics Problem' and Table 1"},{"comment":"The text states that the candidate solution was 'verified and validated using Pontryagin's Principle using the method described in Reference [22]', but no numerical verification results are reported. For a paper whose central claim is a numerically computed solution at unprecedented scale, it is essential to give quantitative verification metrics: e.g., KKT residuals, costate error, or convergence of the solution under mesh refinement. Without such metrics, the reader cannot judge whether the million-point plot is a converged optimum or a numerical artifact.","section":"Verification paragraph in Section 'An Illustrative Nonlinear Nonconvex Astrodynamics Problem'"}],"minor_comments":[{"comment":"The text says the statistics are shown for '21 different values of N', but the table lists only N = 2^10, ..., 2^20, which is 11 values. Please correct the wording or the table.","section":"Table 1"},{"comment":"The run times are normalized but no absolute timings or hardware details are given beyond 'MacBook Air'. Absolute timings would help the reader calibrate the speed claim. Also, the stopping tolerance for GMRES and for the outer SQP iteration is never stated.","section":"Figure 6 and Section 'CONSTRUCTION OF A LOW-MEMORY FAST MILLION-POINT LINEAR SOLVER'"},{"comment":"The approximation (35) is justified by visual inspection of a 10-point matrix (34). This is a heuristic; the paper would be stronger if it reported the spectrum of AP^{-1} for representative ∂X f, or a bound on the condition number of the preconditioned system.","section":"Section 'Development of a Special Birkhoff-Centric Preconditioner'"},{"comment":"There are minor typographical issues, e.g., 'Optimizaiton' in Reference [15], 'Bikrhoff' in the Conclusion, and 'T able' in Table 1. A careful proofreading pass is needed.","section":"Throughout"}],"recommendation":"major_revision","confidential_remarks":"The paper appears to be an AAS conference paper built on a sequence of earlier Birkhoff-theoretic papers by the same group. The gap between the derived scalar fixed-endpoint solver and the demonstrated free-time multi-state orbit transfer is substantial and currently prevents verification of the central claim. The FastPinvY error suggests that the manuscript has not been carefully checked against its own definitions. I would not reject outright because the core ideas (FFT-based FastBV, lower-triangular preconditioner) are plausible and potentially useful, and the missing material is, in principle, fixable. However, the revision must supply the missing KKT-system formulation, a corrected and validated linear solver, and quantitative verification results. Without those, the claims cannot be assessed."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Here's the short version. The paper has a genuinely interesting algorithmic idea—FFT-based Birkhoff matrix-vector products plus a lower-triangular Birkhoff-like preconditioner—and if it works, it would be a real step toward the 2017 million-point challenge. But the central claim is not supported by the text. The linear solver, FastLinSol, is derived only for the fixed-endpoint scalar Newton system [I−B_a ∂_X f]δX = x_a b (Eq. 28). The demonstrated orbit transfer (39) has four states, a control, free final time, and endpoint constraints; its KKT system is never written down. Table 1 reports 'Matrix-free KKT solves' for that problem, but the reader cannot verify what system was solved or how the scalar-diagonal machinery extends to the full Jacobian. More concretely, FastPinvY as printed does not invert the preconditioner P = I−\\tilde{B}_a∂_X f: forward substitution for a lower-triangular P requires accumulating w_j ∂_x f ξ_j, but the algorithm accumulates w_j ∂_x f Υ_j and uses a minus sign. That is an internal inconsistency between Eq. (38) and the algorithm.\n\nWhat the paper does well: the FastBV derivation is clear and mostly self-contained. The observation that the Birkhoff matrix is lower-left dominant and can be approximated by a triangular matrix with Clenshaw-Curtis weights is clever, and the reported constant Krylov iteration count across three orders of magnitude in N is exactly what you would want from a well-preconditioned problem. The authors also honestly scope the novelty: the FFT/DCT trick is standard, and the Birkhoff theory comes from their own earlier papers.\n\nThe soft spots are proportional. The missing KKT system is load-bearing; the million-point claim depends on it. The FastPinvY bug may be a typo, but without code or data we cannot tell. No absolute wall-clock times, no verification metric for the Pontryagin check, and no reproducibility artifacts are provided. The preconditioner is tested on one problem, and the O(1) iteration count is an empirical observation, not a proven bound.\n\nBottom line: this is a paper for the numerical optimal control community. The claim is important enough that it deserves a serious referee rather than a desk reject, and a revised version with the full KKT derivation, a corrected FastPinvY, and code or data could be a solid contribution. As it stands, I would not rely on the million-point result.","headline":"A potentially important solver that, as written, does not demonstrate its headline result: the linear algebra is derived for a simplified fixed-endpoint problem, the actual KKT system is absent, and the printed preconditioner algorithm is inconsistent with its own definition.","tokens_in":20895,"tokens_out":7876,"would_cite":false,"duration_ms":74764,"reading_group":"maybe","serious_thinker":"no","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["49M15","65F08","65T50"],"pacs":[],"model":"deepseek-v4-flash","headline":"This paper presents a fast solver that solves a nonlinear, nonconvex orbit-transfer trajectory optimization over a million grid points on an ordinary laptop, with O(N log N) Birkhoff matrix-vector products and mesh-independent iteration cou","keywords":["trajectory optimization","Birkhoff pseudospectral methods","matrix-free linear algebra","fast Fourier transform","preconditioning","optimal control","orbit transfer","million-point discretization"],"falsifier":"Write down the exact KKT system for the orbit-transfer problem including the free final time and the endpoint constraints, form the preconditioned operator AP^{-1} used by FastLinSol, and count Krylov iterations at N = 1024 and N = 2^20. If the iteration counts grow with N, or if the number of matrix-free KKT solves in the actual run differs from the reported 37, the central mesh-independence claim collapses. A second check: run the same solver on a problem whose dynamics couple adjacent grid points so that ∂X f is not diagonal; if FastPinvY cannot be applied or iterations grow, the method is","tokens_in":19888,"feed_emoji":"🚀","tokens_out":8574,"duration_ms":89699,"temperature":0.7,"pith_summary":"Trajectory optimization normally cannot reach a million grid points: pseudospectral differentiation matrices are dense and ill-conditioned, and off-the-shelf optimizers require explicit Jacobians. This paper claims to break that barrier for a nontrivial, nonlinear, nonconvex orbit-transfer problem by discretizing with a Birkhoff-Chebyshev method whose matrices are well-conditioned and never formed. The key algorithmic step is computing Birkhoff matrix-vector products in O(N log N) time with a discrete cosine transform, and pairing them with a lower-triangular Birkhoff-like preconditioner built from Clenshaw-Curtis weights. On grids from 2^10 to 2^20+1 points, the reported SQP iteration count is constant at 4 and the number of matrix-free KKT solves is constant at 37, yielding a solution on a small laptop in practical time. If correct, this makes dense, high-accuracy pseudospectral grids a realistic tool for optimal control problems that previously forced low-order knotting methods or infeasible memory use.","feed_headline":"Million-point orbit transfer solved fast on a laptop","feed_subtitle":"FFT-based Birkhoff acceleration keeps iteration counts flat as the grid grows, making million-point pseudospectral solves practical.","key_machinery":"The Birkhoff-Chebyshev matrix Ba maps virtual derivative values to state values at Chebyshev-Gauss-Lobatto nodes. Instead of storing Ba, FastBV evaluates Ba V through a DCT-I in O(N log N), using the Chebyshev aliasing formula TN+1 = TN-1 on the grid to keep the mode count consistent. The companion object is the lower-triangular ~Ba defined by Clenshaw-Curtis weights, with diagonal entries half the weight and subdiagonal entries the full weight, which closely approximates Ba and yields the preconditioner P = I - ~Ba ∂X f; P^-1 is applied in O(N) via forward substitution. Together these give an O(N log N) matrix-free linear solver whose Krylov iteration count stays flat as N grows.","core_discovery":"The central discovery is that the Birkhoff-Chebyshev discretization, long known for its O(1) condition numbers, can be made matrix-free at a scale no pseudospectral trajectory solver has reached. The paper constructs algorithms FastBV, FastAX, FastPinvY, and FastLinSol: FastBV applies the dense Birkhoff matrix Ba to a vector by converting to Chebyshev modal coefficients with a DCT-I, using the Chebyshev aliasing formula to absorb the extra mode, for O(N log N) work. FastAX applies the Newton operator I - Ba ∂X f in the same complexity. The preconditioner replaces Ba by a lower-triangular matrix ~Ba whose entries are Clenshaw-Curtis quadrature weights, so P = I - ~Ba ∂X f can be inverted in O","pith_inferences":["The paper derives the fast linear solver for a simplified fixed-endpoint Newton system, then applies it to a free-final-time, endpoint-constrained problem; the actual KKT system is never written. A careful reader should verify that the reported 37 solves are for that full KKT system and not just for the simplified equation.","Because FastAX and FastPinvY rely on ∂X f being diagonal, the O(N log N) guarantee is tied to dynamics whose Jacobian is pointwise diagonal. Problems with path constraints or nonlocal coupling would require block preconditioning, and the O(1) iteration count may degrade.","A natural stress test is a bang-bang control problem with many switches: the near-perfect approximation of Ba by ~Ba relies on smoothness, so hard discontinuities could reveal N-dependent iteration growth that the orbit-transfer example does not.","If the mesh-independence survives those tests, the solver becomes an enabling block for control-oriented orienteering and traveling-salesperson problems that motivated the million-point challenge."],"forward_implications":["A million-point grid is no longer a memory barrier: instead of storing an ~10^12-entry differentiation matrix, the solver stores only O(N) vectors, so high-resolution pseudospectral solutions fit on laptops.","Total solve time is dominated by O(N log N) transforms per Krylov iteration rather than O(N^2) or O(N^3) linear algebra, so the practical cost advantage over direct solvers grows with N.","If iteration counts are genuinely mesh-independent, trajectory optimization cost decouples from grid resolution: refining from 10^3 to 10^6 points does not change the number of expensive optimizer iterations.","The combination of exponential convergence and O(1) conditioning removes the traditional reason for hp/knotting methods, potentially simplifying adaptive mesh refinement in optimal control.","The same matrix-free preconditioning pattern should transfer to other Birkhoff variants and to differential-equation-constrained optimization beyond trajectory problems."],"supporting_citations":[{"why":"Poses the million-point trajectory optimization challenge that this paper answers.","marker":"[1]"},{"why":"Supplies the universal Birkhoff discretization whose condition numbers are O(1), the starting point for the fast solver.","marker":"[8]"},{"why":"Exemplifies the off-the-shelf optimizers that require explicit Jacobian and Hessian data, which the method must avoid.","marker":"[12]"},{"why":"Provides the matrix-free trust-region SQP framework that the new solver builds on for KKT solves.","marker":"[14]"},{"why":"Gives the Hessian and complexity analysis that motivates never forming dense Jacobians or Hessians.","marker":"[15]"},{"why":"Supplies the DCT/FFT and Chebyshev aliasing identities behind the O(N log N) matrix-vector products.","marker":"[18]"},{"why":"Gives the Birkhoff discretization of the orbit-transfer problem and the Pontryagin-based verification used in the demonstration.","marker":"[22]"},{"why":"Defines the Birkhoff basis matrices and Clenshaw-Curtis integrations used to build FastBV and the preconditioner.","marker":"[37]"}],"fun_headline_variants":["Million-point trajectory solves in O(1) iterations","FFT-powered Birkhoff solver reaches million-point scale","Fast million-point trajectory optimization via Birkhoff","Million-gridpoint optimal control with O(1) iteration count"],"cache_read_input_tokens":2688,"weakest_assumption_plain":"The fast solver is derived for a simplified fixed-endpoint version of Newton's equation, and the paper then assumes, without showing the KKT system, that the same O(N log N) matrix-free products and the 37-solve mesh-independence carry over to the full free-final-time, endpoint-constrained orbit-transfer problem.","fun_headline_variants_meta":{"raw":{"variants":["Million-point trajectory solves in O(1) iterations","FFT-powered Birkhoff solver reaches million-point scale","Fast million-point trajectory optimization via Birkhoff","Million-gridpoint optimal control with O(1) iteration count"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000687,"raw_usage":{"total_tokens":2936,"prompt_tokens":713,"completion_tokens":2223,"prompt_tokens_details":{"cached_tokens":256},"prompt_cache_hit_tokens":256,"prompt_cache_miss_tokens":457,"completion_tokens_details":{"reasoning_tokens":2158}},"tokens_in":457,"tokens_out":2223,"duration_ms":23760,"temperature":1.0,"reasoning_tokens":2158,"cache_read_input_tokens":256,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-05T12:07:06.967956+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Write down the exact KKT system for the orbit-transfer problem including the free final time and the endpoint constraints, form the preconditioned operator AP^{-1} used by FastLinSol, and count Krylov iterations at N = 1024 and N = 2^20. If the iteration counts grow with N, or if the number of matrix-free KKT solves in the actual run differs from the reported 37, the central mesh-independence claim collapses. A second check: run the same solver on a problem whose dynamics couple adjacent grid points so that ∂X f is not diagonal; if FastPinvY cannot be applied or iterations grow, the method is","supporting_citations":[{"cited_title":"The Million Point Computational Optimal Con trol Challenge,","cited_arxiv_id":null,"evidence_quote":"Poses the million-point trajectory optimization challenge that this paper answers."},{"cited_title":"A Universal Birkhoff Theory for Fast Traject ory Optimization,","cited_arxiv_id":null,"evidence_quote":"Supplies the universal Birkhoff discretization whose condition numbers are O(1), the starting point for the fast solver."},{"cited_title":"SNOPT: An SQP Al gorithm for Large-Scale Constrained Optimization,","cited_arxiv_id":null,"evidence_quote":"Exemplifies the off-the-shelf optimizers that require explicit Jacobian and Hessian data, which the method must avoid."},{"cited_title":"A Matrix-Free Trust- Region SQP Method for Equality Constrained Optimization,","cited_arxiv_id":null,"evidence_quote":"Provides the matrix-free trust-region SQP framework that the new solver builds on for KKT solves."},{"cited_title":"Hessians in Birkhoff-Theoretic Trajector y Optimizaiton,","cited_arxiv_id":null,"evidence_quote":"Gives the Hessian and complexity analysis that motivates never forming dense Jacobians or Hessians."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Supplies the DCT/FFT and Chebyshev aliasing identities behind the O(N log N) matrix-vector products."},{"cited_title":"Implementations of the Univ ersal Birkhoff Theory for Fast Trajectory Optimization,","cited_arxiv_id":null,"evidence_quote":"Gives the Birkhoff discretization of the orbit-transfer problem and the Pontryagin-based verification used in the demonstration."},{"cited_title":"A Universal Bir khoff Pseudospectral Method for Solv- ing Boundary V alue Problems,","cited_arxiv_id":null,"evidence_quote":"Defines the Birkhoff basis matrices and Clenshaw-Curtis integrations used to build FastBV and the preconditioner."}],"review_version":1}