{"id":"a7431077-cfa1-45d7-88e0-d473413f352c","arxiv_id":"2505.20818","paper_version":1,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":4.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":4,"one_line_summary":"A domain-decomposition subspace neural network method solves linear and nonlinear PDEs with errors down to 1e-13 and lower training cost than PINN, DGM, DRM, and LocELM on 1D/2D benchmarks.","lead":"This paper describes a method that splits a difficult math problem into smaller regions, solves each region with a small neural network, and smoothly stitches the regional answers together. The method reports extremely accurate solutions to test equations while using much less computation than common neural network approaches.","discovery_kind":"extension","skeptic_critique":{"model":"deepseek-v4-flash","headline":"The advertised 1e-11-to-1e-13 accuracies depend on a least-squares solve over basis functions trained only with beta fixed to 1 (Eq. 3); no rank or conditioning analysis supports the global system (13), so the central accuracy claim is not yet established.","rationale":"The reader's weakest assumption is the right one. The central claim is an empirical accuracy/efficiency claim, and the numerical tables are the evidence. The weakest link in the construction is the separated training strategy: basis functions are obtained before coefficients, with no analysis of the quality of the subspace. Because training loss with beta=1 is a single scalar objective, it cannot by itself guarantee that all M basis directions are useful; the later least-squares solve can only be as good as the subspace permits. The omission of condition numbers is particularly salient: for a linear elliptic problem, if the global collocation matrix were well-conditioned and the subspace contained the exact solution, one would expect high accuracy, but the paper does not show either. The ablation tables do provide some independent support: variation over seeds and sampling methods is small, which argues against extreme sensitivity for the tested problems, and hidden-layer ablations show the method is not tuned to a single architecture. But these ablations all keep beta=1 and the same training objective, so they do not test the load-bearing assumption. This is not an internal inconsistency or a demonstrated failure; it is an unverified precondition. A controlled retraining experiment would either confirm the robustness of beta=1 or reveal that the high-accuracy results are an artifact of training along one arbitrary direction. Therefore the appropriate verdict remains conditional, with this check as a required addition. No code is released, so independent reproduction is not currently possible; adding code and the condition-number diagnostic would materially strengthen the claim.","tokens_in":22010,"tokens_out":8468,"duration_ms":94390,"concrete_test":"Settle the concern by rerunning the 1D Helmholtz experiment (Section 5.1.1) with seed 202 and all hyperparameters fixed, but replacing the beta=1 training vector in Eq. (3) with a random fixed vector beta~Uniform[0.5,1.5] (e.g., 5 independent draws), and also reporting cond(A) for the least-squares matrix in (13). If the final L2 error moves outside the 1e-10-1e-12 range or cond(A) exceeds 1e10, the advertised accuracy is not robust to the arbitrary beta=1 training choice.","verdict_should_be":"UNCHANGED","load_bearing_attack":"Section 2 trains the subspace basis by minimizing the scalar PDE loss (3) with beta set to 1, then Section 3 obtains the solution coefficients by least-squares solving (13)-(15). This two-stage construction is the load-bearing path to the headline errors. The training loss constrains only the single all-ones direction of the M-dimensional subspace; the remaining M-1 directions are not explicitly trained, and the stopping rule (4) accepts only a 1e-3 relative residual reduction. Thus all advertised accuracy must be produced by the final least-squares solve, but the paper never reports the numerical rank, smallest singular value, or condition number of the global matrix A in (13), nor whether the basis functions are linearly independent. If A is ill-conditioned or the basis has collapsed directions, the reported 1e-11-1e-13 errors (Tables 1, 6, 21, 25) would be a property of the chosen seed or collocation grid, not of the method. The seed and sampling ablations in Tables 5, 10, 15, and 24 show stability for the tested configurations but do not test the beta=1 choice itself, which is exactly the arbitrary training decision that determines the subspace.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"The manuscript proposes DD-SNN, a domain-decomposition variant of the subspace neural network (SNN) method for linear and nonlinear PDEs. Local SNNs are trained independently on each subdomain by minimizing the PDE residual with the combination coefficients fixed to 1; afterwards, the coefficients are recovered by a least-squares solve of a global system that enforces the PDE, boundary conditions, and C^k interface continuity. Nonlinear problems are handled by Picard and Newton iterations on the coefficients. The paper reports extensive 1D/2D benchmarks (linear and nonlinear Helmholtz, Poisson, parabolic, singular perturbation, viscous Burgers) and claims errors down to 1e-13 with substantially lower training cost than PINNs, DGM, DRM, ELM, and LocELM.","tokens_in":22312,"tokens_out":5537,"duration_ms":59847,"significance":"If the numerical claims hold, the method is a practical contribution: the per-subdomain training is embarrassingly parallel, the coefficient solve is a standard least-squares problem, and the reported errors on smooth manufactured solutions are impressively low. The paper's strengths are the breadth of ablations (subdomain count, subspace dimension, hidden-layer count, sampling method, random seed) and the consistent reporting of errors, epochs, and CPU times. However, the central accuracy claim is supported only by selected benchmark configurations: there is no convergence or conditioning analysis, the interface continuity order used in the experiments is never stated, and the baseline timing comparisons are not controlled. The contribution is promising but needs additional evidence before the headline claims can be accepted.","major_comments":[{"comment":"The load-bearing path to the advertised accuracy is the two-stage construction: basis functions are trained with beta fixed to 1 and stopping tolerance epsilon=1e-3, and then the global coefficient vector is obtained by the least-squares solve (15). The paper provides no analysis or diagnostic showing that this trained subspace is well conditioned for that solve. In particular, the numerical rank, smallest singular value, or condition number of the matrix A in Eq. (13) is never reported. Without such information, the 1e-11 to 1e-13 errors in Tables 1, 6, 21, and 25 may reflect the chosen collocation grid or initialization rather than a property of the method. Please add a conditioning/rank analysis of A for the reported configurations, or an ablation that varies the fixed beta value during basis training and reports the resulting conditioning and errors.","section":"Section 2 (Eq. (3)) and Section 3 (Eqs. (13)-(15))"},{"comment":"The continuity order k_s is a free parameter of the method, but no numerical experiment in Section 5 states which k_s was used. Since the global system (13) is assembled from interface equations of the form (12), the reported accuracy depends directly on this choice, and the omission makes the experiments irreproducible. Please state k_s for each benchmark and justify the choice from the PDE order; the remark following Eq. (5) that a highest derivative k_s requires C^{k_s-1} also needs to be reconciled with the C^k notation used elsewhere.","section":"Section 2 (Eq. (5)) and Section 5"},{"comment":"The comparison baselines do not currently support the claim of 'significantly reducing computational costs.' PINNs, DGM, and DRM are run for a fixed 50,000 epochs without tuning, and CPU times are measured across different implementations (e.g., DEEPXDE for PINNs), so the time comparisons may reflect implementation and hardware rather than algorithmic efficiency. Please provide convergence curves for the baselines, report whether they have plateaued, use equivalent implementations/hardware for all methods, and tune baseline hyperparameters or justify why the chosen settings are representative.","section":"Section 5 (Tables 1, 6, 11, 16, 21, 25)"},{"comment":"The claim that errors 'reach up to 10^-13' overstates the evidence. The 1e-13 accuracy appears only for selected configurations (e.g., Table 2 with 16 subdomains and Table 8 with zero hidden layers), while the singular-perturbation problem in Table 16 saturates around 1e-7 and the nonlinear examples in Tables 21 and 25 are at roughly 1e-9 to 1e-10. The abstract should report the range of errors across benchmarks or state the conditions under which 1e-13 is achieved, and the conclusion should not imply that all tested problems reach that accuracy.","section":"Abstract and Section 5"}],"minor_comments":[{"comment":"The normalization formula appears garbled: the expression '2 bkx−akx · (x − akx) − 1 2 bk t −ak t · (t − ak t ) − 1' is missing parentheses and superscripts, and it is unclear how the temporal coordinate is treated in steady-state problems where no time direction exists.","section":"Eq. (2a)"},{"comment":"The text says DD-SNN achieves its result with 'an average of 135 epochs,' while Table 11 reports 98 epochs for the same configuration; please correct the inconsistency.","section":"Section 5.1.3 (Table 11)"},{"comment":"For the two-dimensional Poisson equation with 16 subdomains and uniform sampling, Table 6 reports an L_infinity error of 6.61e-11 while Table 9 reports 2.02e-10 for what appears to be the same configuration; please verify which value is correct.","section":"Tables 6 and 9"},{"comment":"The text states that the PINN L2 error is 9.96e-3, but Table 11 reports 2.29e-2 for PINN; please correct the text or the table.","section":"Section 5.1.3 (Table 11)"},{"comment":"Reference [19] is cited for two different works in the introduction: once for 'compensated deep Ritz methods' and once for 'randomized neural network with Petrov-Galerkin methods'; please renumber the references so each citation points to the intended work.","section":"References in Section 1"}],"recommendation":"major_revision","confidential_remarks":"The paper would benefit from an editor check on the fairness of the DGM/DRM baselines, since these methods are usually deployed with problem-specific tuning and the fixed 50k-epoch runs may not be representative. The missing k_s specification is a straightforward but important reproducibility fix. The core idea is worth a revision rather than rejection."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"What you should know first: this is a solid engineering paper, not a theory paper. The method extends the SNN and LocELM lines by training local subspace neural-network bases per subdomain with only the PDE residual (beta fixed to 1), then solving a global least-squares problem for the combination coefficients with C^k interface constraints. Nonlinear equations are handled by Picard or Newton iteration on the coefficients alone. The numerical study is unusually thorough: for six problems (linear and nonlinear, steady and time-dependent) they sweep subspace dimension, subdomain count, hidden-layer depth, sampling scheme, and random seed, and report errors down to 1e-13.\n\nThe genuinely new piece is the specific combination—subspace basis training plus domain decomposition plus C^k interface conditions plus coefficient-only nonlinear iterations. That combination does not appear in the cited SNN [8] or LocELM [10] papers, and the ablation tables show that increasing subdomain count generally reduces both error and training epochs, up to a point. The method is plausible and the core idea is sound.\n\nThe soft spots are in the evidence, not the concept. The load-bearing accuracy claim depends on a least-squares solve over bases trained only along the beta=1 direction, but the paper reports no condition number, rank, or smallest singular value for the global matrix in (13). We do not know whether the 1e-11 errors are a property of the method or of the chosen seed and collocation grid. The seed and sampling ablations mitigate this somewhat—six seeds all give roughly 1e-11 on the Poisson problem—but they do not test the beta=1 training choice itself. Missing too: any convergence analysis, the actual C^k order used in the experiments, and repeated-run statistics. There is also an internal inconsistency in the 2D Poisson results: Table 6 reports L_infinity = 6.61e-11 for the baseline, while Table 9 lists 2.02e-10 for the same configuration. The comparison baselines (PINN, DGM, DRM) all get fixed 50,000 epochs with no tuning, so the 'significantly reduces computational cost' claim is fair in this controlled setting but not against well-tuned competitors. No code is released, so the headline numbers are not independently checkable.\n\nWorth sending to a serious referee? Yes. It is a legitimate extension with strong numerical evidence, and the missing pieces—conditioning analysis, convergence insights, clarified ks, consistent tables, and code—are exactly what a referee should request. My recommendation: engage with it, but push for a revision that addresses the conditioning of the global system before the accuracy claims are taken at face value.","headline":"A useful engineering combination of domain decomposition, subspace neural bases, and coefficient-only nonlinear iterations, with extensive but not yet fully supported accuracy claims.","tokens_in":22783,"tokens_out":3051,"would_cite":false,"duration_ms":28964,"reading_group":"maybe","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["65N55","65N35","68T07"],"pacs":[],"model":"deepseek-v4-flash","headline":"A neural PDE solver that trains local basis functions first, then fits their combination coefficients under interface smoothness constraints, reaches errors near 1e-13 while using far fewer epochs than end-to-end physics-informed networks.","keywords":["domain decomposition","subspace neural network","partial differential equations","C^k continuity","least-squares coefficient solve","Picard iteration","Newton iteration","neural network PDE solvers"],"falsifier":"On the paper's 1D Helmholtz benchmark, train the local bases with $\\beta=1$ on one set of collocation points, determine the coefficients on a second set, and evaluate the assembled solution on a dense independent grid, recording the $L^2$ error and the condition number of the least-squares matrix; if the independent error is far above the reported $10^{-11}$ level, or if the matrix condition number grows sharply with subspace dimension, the separated-training premise is not enough to support the claimed accuracy.","tokens_in":21829,"feed_emoji":"🧮","tokens_out":11150,"duration_ms":111464,"temperature":0.7,"pith_summary":"This paper proposes a two-stage method for PDEs: partition the domain into non-overlapping subdomains, train a small neural network in each subdomain to produce a set of basis functions by minimizing only the PDE residual (with combination coefficients fixed to 1), then freeze the bases and solve a global least-squares problem for the coefficients, imposing $C^k$ continuity across interfaces. The paper's claim is that this scheme solves linear and nonlinear benchmark problems — Helmholtz, Poisson, parabolic, singular-perturbation, and viscous Burgers — with errors of order $10^{-10}$ to $10^{-13}$ while using hundreds rather than tens of thousands of training epochs. If true, this matters because it indicates that global end-to-end training is not what limits neural PDE accuracy; a local residual fit plus a constrained coefficient solve can carry most of the work. The paper also claims that nonlinear equations can be handled by Picard or Newton iteration on the coefficients alone, with the basis fixed.","feed_headline":"Basis-first neural solver hits PDE errors near 1e-13","feed_subtitle":"Domain decomposition plus a least-squares coefficient fit beats PINNs on linear and nonlinear benchmarks","key_machinery":"The key object is the trained local basis set. Each subdomain network ends in a subspace layer of dimension $M$ that outputs basis functions $\\phi_1,\\ldots,\\phi_M$, and the local solution is a linear combination $\\hat u = \\beta\\cdot\\Phi$. The mechanism is the two-stage split: first train the network weights by minimizing the PDE residual with $\\beta$ fixed to 1; then freeze the basis and solve an overdetermined least-squares system built from rows that enforce the equation at collocation points, the boundary conditions, and $C^k$ continuity (matching function values and derivatives up to order $k_s$) on shared interfaces. That global system, not the neural training, produces the smooth assembled solution and provides the linearized updates used in Picard and Newton iterations for nonlinear equations.","core_discovery":"The discovery the paper is trying to establish is that separating basis generation from solution fitting makes high-accuracy neural PDE solving cheap. In each subdomain, a subspace neural network outputs $M$ basis functions and represents the local solution as $\\hat u = \\sum_{j=1}^M \\beta_j \\phi_j$; training updates only the network weights, with $\\beta$ fixed to 1, and then the assembled least-squares system determines $\\beta$ from three sources: interior equation residuals, boundary conditions, and $C^k$ interface conditions between adjacent subdomains. On the reported benchmarks this two-step fit reaches $L^2$ errors of $4.98\\times10^{-11}$ (1D Helmholtz), $4.94\\times10^{-12}$ (2D Poisson), $2.33\\times10^{-11}$ (parabolic), $2.45\\times10^{-7}$ (singular-perturbed boundary layer), and $9.49\\times10^{-10}$ (nonlinear Helmholtz), with a best 1D Helmholtz error of $9.49\\times10^{-13}$ at 16 subdomains, and the Newton version on viscous Burgers reaches $4.19\\times10^{-10}$. The paper interprets these numbers as showing that domain decomposition and subspace representation, not network depth or long training, are the decisive accuracy ingredients.","pith_inferences":["A likely reading, not stated in the paper, is that most of the accuracy comes from the overdetermined least-squares fit and the interface constraints rather than from the trained hidden layers; the paper's own tables show that zero hidden layers is sometimes the most accurate configuration, which makes a clean testable comparison between trained bases and untrained random-feature bases.","The $C^k$ interface condition scales with equation order: a fourth-order PDE would require $C^3$ continuity across interfaces, adding derivative rows to the global system and potentially eroding the efficiency advantage reported here.","The reported errors are all against smooth, manufactured exact solutions; a natural stress test is a nonsmooth or shock-forming solution where a smooth residual-trained basis may require many subdomains or fail to represent the discontinuity.","The two-stage structure invites an adaptive loop, not present in the paper: after the coefficient solve, measure the residual locally and refine subdomains or enlarge the subspace where the residual is largest."],"forward_implications":["For a fixed trained basis, changing the source term or boundary data in a linear problem reduces to re-solving the least-squares system for the coefficients, so parametric sweeps avoid retraining the networks.","Because each subdomain trains independently, the method supports straightforward parallel training, and increasing resolution can be done by adding subdomains rather than deepening the network.","For nonlinear problems, Picard and Newton iterations act on the small coefficient vector with the basis frozen, so per-iteration cost is governed by the least-squares solve rather than by network backpropagation.","The experiments show accuracy improving as the number of subdomains increases up to a limit set by the number of collocation points per subdomain, so domain count and subspace dimension are practical accuracy controls."],"supporting_citations":[{"why":"provides the subspace neural network architecture whose basis-function representation DD-SNN extends to subdomains.","marker":"[8]"},{"why":"supplies the local extreme learning machine baseline and the nonlinear least-squares strategy that DD-SNN replaces.","marker":"[10]"},{"why":"introduces the physics-informed neural network baseline that the paper compares against.","marker":"[5]"},{"why":"gives the non-gradient least-squares iteration approach that the paper's coefficient-update schemes improve on.","marker":"[11]"},{"why":"introduces the earlier domain-decomposition physics-informed network framework used as motivation for interface continuity.","marker":"[27]"}],"fun_headline_variants":["Subspace net split beats PINNs, hits 1e-13","Domain-split basis nets edge out PINNs at 1e-13","Basis-first neural net nails PDEs to 1e-13","Neural subspace solver: 1e-13 errors, fewer training costs"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The load-bearing premise is that basis functions trained with $\\beta=1$ against the bare PDE residual form a subspace rich enough and stable enough for the later least-squares coefficient solve to deliver the advertised $10^{-10}$ to $10^{-13}$ accuracy; the paper provides no analysis connecting that training objective to the final fit.","fun_headline_variants_meta":{"raw":{"variants":["Subspace net split beats PINNs, hits 1e-13","Domain-split basis nets edge out PINNs at 1e-13","Basis-first neural net nails PDEs to 1e-13","Neural subspace solver: 1e-13 errors, fewer training costs"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000653,"raw_usage":{"total_tokens":2986,"prompt_tokens":931,"completion_tokens":2055,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":547,"completion_tokens_details":{"reasoning_tokens":1974}},"tokens_in":547,"tokens_out":2055,"duration_ms":14899,"temperature":1.0,"reasoning_tokens":1974,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-07T13:46:57.462544+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"On the paper's 1D Helmholtz benchmark, train the local bases with $\\beta=1$ on one set of collocation points, determine the coefficients on a second set, and evaluate the assembled solution on a dense independent grid, recording the $L^2$ error and the condition number of the least-squares matrix; if the independent error is far above the reported $10^{-11}$ level, or if the matrix condition number grows sharply with subspace dimension, the separated-training premise is not enough to support the claimed accuracy.","supporting_citations":[{"cited_title":"Subspace method based on neural networks for solving the partial differential equation","cited_arxiv_id":"2404.08223","evidence_quote":"provides the subspace neural network architecture whose basis-function representation DD-SNN extends to subdomains."},{"cited_title":"Raissi, P","cited_arxiv_id":null,"evidence_quote":"introduces the physics-informed neural network baseline that the paper compares against."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"gives the non-gradient least-squares iteration approach that the paper's coefficient-update schemes improve on."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"introduces the earlier domain-decomposition physics-informed network framework used as motivation for interface continuity."}],"review_version":1}