{"id":"08fcfcc3-bf0d-434f-908d-8e6d9f40fbb7","arxiv_id":"1909.00049","paper_version":1,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":7.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":7,"one_line_summary":"A panel-based high-order adaptive boundary integral scheme solves 2D Stokes flow in complex nonsmooth domains to user-specified precision with near-linear cost.","lead":"An adaptive boundary integral method is built for slow viscous (Stokes) flows in complicated 2D geometries with corners, using special quadrature rules for nearly singular integrals. It solves a 378-corner vascular network to 9-digit accuracy, which matters for microfluidic, suspension, and blood flow simulations that need fast, automatic, high-precision solvers.","discovery_kind":"new_method","skeptic_critique":{"model":"deepseek-v4-flash","headline":"Close-evaluation Vandermonde blocks are the unverified linchpin: backward stability alone does not imply accurate near-corner integrals, and nonsmooth accuracy claims rest on self-convergence.","rationale":"The reader's weakest assumption points to the same mechanism, and I agree: the evidence in the paper is strong but internal. The smooth-domain tests give some independent confidence, and the FMM/timing data are plausible; no mathematical contradiction inside the core derivation was found. However, the central nonsmooth accuracy claim rests on the correctness of these ill-conditioned Vandermonde blocks plus the geometric grading. A failure there would not be detected by the paper's self-convergence check because the comparison solution uses the same blocks. This justifies keeping the verdict conditional rather than accepting the headline as established. The O(N^(3/2)) close-correction scaling is a second, independently documented concern about the linear-scaling label, but it is secondary to accuracy and does not change the recommended verdict.","tokens_in":23345,"tokens_out":18447,"duration_ms":180252,"concrete_test":"Implement a standalone version of the close-evaluation block (Sections 3.3.1-3.3.4) for a source panel adjacent to the reentrant corner used in Example 2 or 4, with target nodes on the neighboring panel across the corner. Use a known exact scalar density (e.g., tau(y)=1 and tau(y)=y^3) and compare the resulting IL, IC, IH, IS values against the same contour integrals evaluated with high-precision adaptive quadrature (50-digit arithmetic). Repeat for the m values used in the paper (m=2p, p=3..13). If any relative error at the requested tolerance (say 1e-9) exceeds 1e-9, the Vandermonde/backward-stability assumption fails exactly in the regime the paper relies on.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The decisive link in the argument is the close/self-evaluation quadrature of Section 3.3, used to fill the near-diagonal blocks of the Nyström matrix. For the Cauchy, logarithmic, Hadamard and supersingular potentials, the derived scalar density is represented by a degree m-1 polynomial in the complex variable y, and the coefficients are obtained by solving the Vandermonde system (29) (via the adjoint solve V^T X = P^T in Section 3.4). The system is exponentially ill-conditioned in m; backward stability of mldivide only controls a residual of the solved linear system, not the error in the resulting integrals P V^{-1} L when V is ill-conditioned and the monomial coefficients are large. The paper's assurance that m<50 is safe is a heuristic taken from [27, App. A], and the paper supplies no independent check on a corner-adjacent panel pair. This matters because all nonsmooth convergence evidence in Section 5 is self-referential: the exact solution is not known for the corner/vascular examples, so a systematic error in these blocks would be invisible in the reported convergence to the finest grid solution but would invalidate the claimed 1e-9 accuracy.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"The manuscript presents an adaptive, panel-based Nyström boundary integral scheme for the Stokes Dirichlet problem in two dimensions on complex, multiply connected, nonsmooth domains. The Stokes velocity potentials are rewritten in terms of four scalar complex contour integrals (logarithmic, Cauchy, Hadamard, and supersingular), for which close- and self-evaluation quadratures are constructed by polynomial interpolation of the derived scalar densities and exact monomial integration. Self- and close-touching targets are treated by a unified rule with branch-cut rotation. A three-stage adaptive panel refinement algorithm chooses panel degree p and panel lengths from a user tolerance, with geometric corner grading and closeness refinement. The paper reports superalgebraic convergence on a smooth starfish example, root-exponential convergence for the shuriken and multi-polygon examples, and solves a 378-corner vascular network to a reported relative L2 error of 1e-9 with 356,580 unknowns, with GMRES iteration counts and timing breakdowns. Pressure and traction formulas are derived in an appendix and validated on the smooth example.","tokens_in":23597,"tokens_out":10920,"duration_ms":103223,"significance":"If the claims are sustained, this is a useful practical contribution: it directly handles Stokes layer potentials in physical variables, unifies on-surface and near-surface evaluation, provides a user-tolerance-driven adaptive discretization for nonsmooth multiply connected geometries, and includes pressure and traction variants. The strength of the paper is its combination of detailed quadrature derivation, a concrete adaptive algorithm, and extensive numerical experiments with convergence tables and timings for five distinct geometries. The weak points are that the corner-accuracy claims are measured only against a finer run of the same scheme and that the close-evaluation quadrature rests on an ill-conditioned Vandermonde solve whose corner-adjacent behavior is not independently validated; both need to be addressed before the quantitative accuracy claims can be taken at face value.","major_comments":[{"comment":"All reported errors in the nonsmooth examples are measured against the finest-grid solution of the same scheme, as the authors explicitly state ('the exact solution is not known analytically; therefore, we use the finest grid solution as the reference solution'). This makes the claimed 1e-9 relative L2 error in the 378-corner example an internal self-convergence measure, not a verification of true accuracy: a systematic quadrature error that persists across all refinements, for instance in the close-evaluation blocks of Section 3.3, would be invisible. Please add at least one nonsmooth test with an independent reference, such as a wedge or polygon problem with a known exact Stokes solution containing the corner-singular terms, or a comparison against an independent high-accuracy method (for example, RCIP or an over-resolved finite-element benchmark), and report both the self-convergence error and the error against that independent reference.","section":"Section 5, Table 1 and Examples 2–4"},{"comment":"The close-evaluation quadrature rests on solving the Vandermonde system V a = tau_tilde, and the matrix block A is formed via the adjoint solve V^T X = P^T. The paper justifies this with the statement, following [27, App. A], that backward stability of MATLAB mldivide suffices for m<50. Backward stability bounds the residual of the solved linear system, not the forward error in A = X^T L when V is exponentially ill-conditioned and the polynomial coefficients are large. Since all near-corner evidence in Section 5 is self-referential, a systematic error in exactly these blocks would be invisible. Please provide a targeted validation: for a panel touching a reentrant corner and a neighboring target panel (or on-surface target near the corner), compare the Cauchy, logarithmic, Hadamard, and supersingular outputs of A against a high-precision direct quadrature (for example, 64-digit arithmetic), and report sensitivity to m = p, 2p, 3p. If the m<50 heuristic does not hold for corner-adjacent panels, the adaptive corner refinement must be modified.","section":"Section 3.3.1, Eq. (29) and Section 3.4"},{"comment":"The paper claims a linear-scaling scheme and states in the Conclusions that 'CPU time grows linearly with problem size', but its own asymptotic estimate in Section 5 states that the close-evaluation matrix-vector multiplication time grows as O(N^{3/2}), and Table 1 shows the FMM fraction decreasing from 78.5% to 66.2% as epsilon decreases. Thus the total per-iteration cost is not asymptotically linear, and the claim that the scheme is linear-scaling is not supported beyond the tested tolerance range. Please either implement an O(N) close-correction application (for example, by exploiting low-rank or translation-invariant structure) or revise the title and conclusion claims to 'near-linear for moderate tolerances', with the O(N^{3/2}) close-evaluation bottleneck stated in the abstract and conclusions.","section":"Title, Section 5 (Table 1 discussion), and Conclusions"}],"minor_comments":[{"comment":"In the displayed formula for q_k, the intermediate expression contains log(1+x), while the derivation and the subsequent branch-cut form use -1-x; please correct the notation to avoid ambiguity.","section":"Section 3.3.2, Eq. (40)"},{"comment":"The claim that flipping the sign of phi is robust for corners of arbitrary angles is not supported by an explicit experiment; a small parameter study sweeping the interior corner angle between 0 and 2pi would strengthen the paper.","section":"Section 3.3.1, Remark 3"},{"comment":"After presenting the adjoint construction for the Cauchy kernel, please state explicitly that the same construction is used for the logarithmic, Hadamard, and supersingular blocks, and define n (the number of targets) consistently in all cases.","section":"Section 3.4"},{"comment":"Please add a legend to Figure 8 and state whether the horizontal axis is the number of nodes N or the number of unknowns 2N, since the text refers to 'required number of unknowns' without specifying the convention.","section":"Figure 8"},{"comment":"The statement 'approximately 800 degrees of freedom per corner' should be checked against the conventions used in Table 1 and Example 4; the ratio 356,580/378 = 943 uses 2N, whereas the text elsewhere sometimes counts N or 2N without consistent definition.","section":"Section 5, Example 3"}],"recommendation":"major_revision","confidential_remarks":"The paper is within the scope of the journal and the numerical machinery is promising, but the quantitative accuracy claims currently rely on self-convergence and on an unvalidated ill-conditioned Vandermonde step. I recommend requiring an independent nonsmooth benchmark and a direct validation of the close-evaluation blocks near corners before publication. No concerns about citation practice or novelty; the self-citations to [25,27] are appropriate and the relationship to Ojala-Tornberg [26] is clearly stated."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Short version: this is a genuinely useful methods paper, and the main algorithm probably works. It is the first Nyström panel scheme I know that handles all the Stokes layer potentials (velocity, pressure, traction) on nonsmooth multiply connected 2D domains with adaptive panel refinement. The core quadrature is a direct extension of Helsing–Ojala and the authors' own smooth-curve scheme, but the extension to arbitrary panels with corners, plus the adaptive refinement pipeline, is real new work.\n\nWhat's good: the paper is clearly written, the quadrature derivations are careful, and the numerical evidence is substantial. The 378-corner vascular network example reaching 1e-9 with 356,580 unknowns is impressive. The authors are honest that the refinement rules are ad hoc; they do not oversell the theory, and they cite the hp-BEM corner complexity work appropriately. The close-evaluation formulas for pressure and traction in the appendix are a useful addition.\n\nSoft spots, in order of importance. First, no code or data is released. For a numerical scheme of this intricacy, independent reproduction is the only way to trust the corner claims, and the paper as it stands cannot be checked. Second, the corner convergence evidence is self-referential: exact solutions are unknown, and the 1e-9 claim is relative to the finest-grid solution. The stress-test note about the Vandermonde close-evaluation blocks is the sharpest version of this worry. The paper relies on backward stability of MATLAB's mldivide for an exponentially ill-conditioned Vandermonde system with m<50; that is plausible and follows [27], but there is no standalone test on a corner-adjacent panel pair against a known integral. A systematic error in those blocks would be invisible in the convergence plots. I do not think this is fatal, but it should be addressed, e.g. with a manufactured wedge solution or a direct reference integral on a corner panel. Third, 'linear-scaling' in the title is stronger than the evidence: Table 1 shows close-correction time growing like O(N^{3/2}) in the tight-tolerance regime, and the FMM dominates only at large N. The conclusion's 'CPU time grows linearly' should be reworded. Fourth, the many tuning constants (C, alpha, lambda, beta, c', phi) are empirically chosen; a sensitivity study would help, though this is minor.\n\nBottom line: I would send this to a serious referee. The core scheme is sound and likely to be used; the paper deserves publication after code release and a non-self-referential corner test.","headline":"A genuinely useful and probably sound adaptive Nyström BIE scheme for nonsmooth 2D Stokes flow, with the main caveat that the corner-accuracy claims rest on self-convergence and need an independent check.","tokens_in":24151,"tokens_out":2915,"would_cite":true,"duration_ms":28046,"reading_group":"yes","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["65N38","76D07","65R20"],"pacs":[],"model":"deepseek-v4-flash","headline":"This paper develops an adaptive boundary integral method that solves two-dimensional Stokes flow in complex, corner-laden geometries to user-set accuracy, demonstrated at 1e-9 relative error in a 378-corner vascular network.","keywords":["Stokes flow","boundary integral equations","Nyström method","nearly singular quadrature","close evaluation","adaptive mesh refinement","corner singularities","fast multipole method"],"falsifier":"Take a Stokes Dirichlet problem with a known exact solution on a domain containing a reentrant corner with an unusual angle, disable the corner-grading rule so panels near the corner are not geometrically refined, and then measure the relative L2 error as the polynomial degree increases: if the error stops decreasing once the panel density develops a non-polynomial singularity, the claim that adaptive grading reduces every panel to a polynomial regime is shown to carry the accuracy statement.","tokens_in":23094,"feed_emoji":"🌊","tokens_out":7523,"duration_ms":71227,"temperature":0.7,"pith_summary":"The paper builds a complete numerical pipeline for low-Reynolds-number flow around and through complicated two-dimensional objects: boundary panels are chosen automatically, every boundary integral is evaluated by high-order quadrature even when the evaluation point sits on or arbitrarily near the boundary, and the resulting linear system is solved iteratively. The central claim is that a user can specify an error tolerance and the scheme will meet it without hand-tuning the mesh, with convergence faster than any power law on smooth boundaries and root-exponential convergence (error of order $e^{{-c sqrt N}}$) when corners are present. If true, this makes boundary integral methods practical for microfluidic channels, porous media, and blood-vessel networks, where geometry is complex and surfaces nearly touch. The demonstration case is a 378-corner vascular network solved to relative L2 error 1e-9 with fewer than 200K boundary points.","feed_headline":"Stokes solver hits 9 digits in a 378-corner vessel network","feed_subtitle":"Adaptive panel rules evaluate nearly singular Stokes integrals automatically, hitting user-set tolerances in nonsmooth domains.","key_machinery":"The load-bearing object is the close-evaluation panel quadrature for the four contour integrals I_L, I_C, I_H, and I_S (logarithmic, Cauchy, Hadamard, and supersingular) into which the Stokes single- and double-layer potentials are decomposed. For a target near or on a panel, the derived scalar density is first upsampled to twice the panel degree using fine Legendre nodes, then approximated by a polynomial in the complex variable y via a backward-stable solve of a Vandermonde system; the monomial integrals follow from a recurrence whose base case is the Cauchy integral of a constant, evaluated as a logarithm with a rotated branch cut so the cut is pushed behind the panel. This one construction supplies accurate matrix entries for self-interaction, nearly singular interactions between nearby panels, and off-surface evaluation, and it is what makes the adaptive refinement scheme trustworthy.","core_discovery":"The paper's claim is that all Stokes layer potentials needed for a Dirichlet velocity problem—single layer, double layer, and the associated pressure and traction fields—can be evaluated to high accuracy at targets on or arbitrarily close to the boundary by expressing them as combinations of four complex contour integrals (logarithmic, Cauchy, Hadamard, and supersingular) and applying one unified panel quadrature to each. On each panel the derived scalar density is approximated by a polynomial in the complex coordinate y; the monomial integrals are computed through a recurrence seeded by a logarithm whose branch cut is rotated to lie behind the panel, which removes the need to distinguish close evaluation from self evaluation. Around corners the panels are geometrically graded, and panels are also refined when boundaries approach one another, so the discretization adapts automatically to a user-set tolerance. The paper verifies in a 378-corner model vascular network that the relative L2 velocity error reaches 1e-9 with 356,580 degrees of freedom, with GMRES converging steadily and the fast multipole method dominating the CPU time.","pith_inferences":["The same contour-integral reduction may carry over to other elliptic kernels with Cauchy-type and supersingular contour integral analogues, such as Helmholtz, elasticity, or biharmonic problems, since the quadrature machinery is built on the four integrals rather than on Stokes specifics.","The corner treatment is geometric rather than analytic; introducing a per-corner singular basis or a compression step for corner blocks would likely lower the degrees of freedom per corner and is a natural next step the paper itself flags.","A rigorous a posteriori error estimate for the adaptive refinement criteria would let the scheme quote certified tolerances rather than empirically matched ones; the paper demonstrates tolerance matching numerically but does not prove it.","In time-dependent particulate flows, one could keep the adaptively paneled fixed channel geometry unchanged while moving particles are handled by global close-evaluation quadrature; whether compressed operator representations can be updated as panels change is an open question the paper poses."],"forward_implications":["A user can hand a complex 2D geometry, including corners and near-touching boundaries, to the solver with a tolerance and receive a solution that meets it; no hand-tuned panel distribution is needed.","For smooth boundaries the error converges faster than any power of the number of unknowns, and with corners it converges like exp(-c sqrt N), so a few hundred thousand unknowns can reach 1e-9 even with hundreds of corners.","Because layer potential evaluation is accurate arbitrarily close to the boundary, the same machinery supports computing velocity, pressure, and traction fields, including force calculations on closely interacting bodies.","For two smooth surfaces separated by a distance d, the adaptive criterion costs about O(log(1/d)) unknowns where a uniformly resolved grid costs O(1/sqrt d), making dense suspensions of nearly touching particles feasible.","The near-linear CPU scaling, dominated by the fast multipole application, means the method can serve as the stationary-geometry component of an evolving particulate flow simulation."],"supporting_citations":[{"why":"Supplies the close-evaluation panel quadrature strategy for Laplace potentials that this work extends to Stokes potentials.","marker":"[27]"},{"why":"Provides the prior global close-evaluation scheme for smooth closed curves whose limitations motivate the adaptive panel approach.","marker":"[25]"},{"why":"Gives the theoretical root-exponential convergence rate for hp-boundary element methods on polygons that the numerical results are compared against.","marker":"[37]"},{"why":"Analyzes the power-law behavior of Stokes densities near corners, which informs the choice of corner grading exponents.","marker":"[39]"},{"why":"Provides the fast multipole method for the biharmonic kernel used to apply the dense operator in near-linear time.","marker":"[45]"},{"why":"States the Bernstein-ellipse error bound for Legendre quadrature used to set the global polynomial degree from the tolerance.","marker":"[47]"},{"why":"Documents the exponential ill-conditioning of Vandermonde matrices, supporting the reliance on backward-stable direct solves in the close-evaluation quadrature.","marker":"[49]"}],"fun_headline_variants":["Adaptive Stokes scheme hits 9 digits in 378-corner network","Fast Stokes solver adapts to corners, nailing 9 digits","Linear-scaling Stokes scheme solves complex 2D corner geometry","Adaptive high-order Stokes solver handles 378 corners at 9 digits","Auto-adaptive Stokes quadrature solves nonsmooth 2D flows"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The entire method assumes that on every panel the unknown density behaves like a low-degree polynomial in the panel's curved coordinate, and that the numerical routine used to fit that polynomial stays accurate even though the fitting equations are extremely ill-conditioned.","fun_headline_variants_meta":{"raw":{"variants":["Adaptive Stokes scheme hits 9 digits in 378-corner network","Fast Stokes solver adapts to corners, nailing 9 digits","Linear-scaling Stokes scheme solves complex 2D corner geometry","Adaptive high-order Stokes solver handles 378 corners at 9 digits","Auto-adaptive Stokes quadrature solves nonsmooth 2D flows"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.001338,"raw_usage":{"total_tokens":5441,"prompt_tokens":948,"completion_tokens":4493,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":564,"completion_tokens_details":{"reasoning_tokens":4401}},"tokens_in":564,"tokens_out":4493,"duration_ms":28998,"temperature":1.0,"reasoning_tokens":4401,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-14T10:03:58.677477+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Take a Stokes Dirichlet problem with a known exact solution on a domain containing a reentrant corner with an unusual angle, disable the corner-grading rule so panels near the corner are not geometrically refined, and then measure the relative L2 error as the polynomial degree increases: if the error stops decreasing once the panel density develops a non-polynomial singularity, the claim that adaptive grading reduces every panel to a polynomial regime is shown to carry the accuracy statement.","supporting_citations":[{"cited_title":"Helsing and R","cited_arxiv_id":null,"evidence_quote":"Supplies the close-evaluation panel quadrature strategy for Laplace potentials that this work extends to Stokes potentials."},{"cited_title":"Spectrally accurate quadratures for evaluation of layer potentials close to the boundary for the 2D Stokes and Laplace equations","cited_arxiv_id":null,"evidence_quote":"Provides the prior global close-evaluation scheme for smooth closed curves whose limitations motivate the adaptive panel approach."},{"cited_title":"The hp-version of the boundary element method on polygons","cited_arxiv_id":null,"evidence_quote":"Gives the theoretical root-exponential convergence rate for hp-boundary element methods on polygons that the numerical results are compared against."},{"cited_title":"On the solution of Stokes equation on regions with corners","cited_arxiv_id":"1711.04072","evidence_quote":"Analyzes the power-law behavior of Stokes densities near corners, which informs the choice of corner grading exponents."},{"cited_title":"On the numerical solution of the biharmonic equation in the plane","cited_arxiv_id":null,"evidence_quote":"Provides the fast multipole method for the biharmonic kernel used to apply the dense operator in near-linear time."},{"cited_title":"Approximation theory and approximation practice , volume 128","cited_arxiv_id":null,"evidence_quote":"States the Bernstein-ellipse error bound for Legendre quadrature used to set the global polynomial degree from the tolerance."},{"cited_title":"How bad are Vandermonde matrices? SIAM J","cited_arxiv_id":null,"evidence_quote":"Documents the exponential ill-conditioning of Vandermonde matrices, supporting the reliance on backward-stable direct solves in the close-evaluation quadrature."}],"review_version":1}