Abstract
<title>Abstract</title> <p> We present a robust and adaptive meshless collocation method that uses strong-form radial basis functions (RBFs) to solve nonlinear elliptic boundary value problems on irregular two-dimensional domains. Our approach uses locally scaled polyharmonic-spline RBF stencils with quadratic poly-augmentation, sparse analytic Jacobians, damped Newton steps, and an off-node residual estimator that allows for automatic refinement of the point cloud. By normalizing local coordinates, we address fill-distance scaling in the interpolation matrices and reduce the method's sensitivity to the shape parameters of smooth kernels. For uniform elliptic and strongly monotone quasilinear operators, our method shows conditional a posteriori reliability, based on a computable discrete-stability result. The estimator uses weighted internal residuals and boundary defects, working within a solve-estimate-mark-refine loop that supports adding new centers based on residual ranking and distance regularization. We tested the method on three nonlinear benchmarks: a smooth solution on a five-lobed star, a corner-singular solution on an L-shaped domain, and a strongly localized layer on a multiply connected flower-shaped domain. After four iterations, the nonlinear residual dropped below 1.7 x 10 <sup>-1</sup> in 24 adaptive and quasi-uniform runs. With the same number of nodes, adaptive refinement reduced the relative <italic>L</italic> <sup>2</sup> error by 60.4% in the re-entrant domain and by 67.3% in the localized-layer benchmark. The H¹ error in the localized-layer case dropped by 71.8%. The global estimator had a Pearson correlation of 0.997 with the relative H1 error on a log-log scale. Local scaling increased the median stencil condition numbers by a factor of 12.4 to 19.0. This framework is flexible for different geometries, reproducible, and offers a high-accuracy baseline for nonlinear diffusion, reaction-diffusion, and similar elliptic models where traditional mesh generation is difficult. </p>