A new meshless method for solving nonlinear boundary value problems, based on the local boundary integral equation (LBIE) method and the moving least squares approximation, is proposed in the present paper. The total formulation and a rate formulation are developed for the implementation of the present method. The present method does not need domain and boundary elements to deal with the volume and boundary integrals, which will cause some difficulties for the conventional boundary element method (BEM) or the field/boundary element method (FBEM), as the volume integrals are inevitable in dealing with nonlinear boundary value problems. This is the same for the element free Galerkin (EFG) method which also needs element-like cells in the entire domain to evaluate volume integrals. The "companion fundamental solution" introduced in Zhu, Zhang and Atluri (1998) is used so that no derivatives of the shape functions are needed to construct the stiffness matrix for the interior nodes, as well as for those nodes with no parts of their local boundaries coinciding with the global boundary of the domain of the problem, where essential boundary conditions are specified. It is shown that the satisfaction of the essential as well as natural boundary conditions is quite simple, and algorithmically very efficient, in the present nonlinear LBIE approach. Numerical examples are presented for several problems, for which exact solutions are available. The present method converges fast to the final solution with reasonably accurate results for both the unknown variable and its derivatives. No post processing procedure is required to compute the derivatives of the unknown variable (as in the conventional FBEM), since the solution from the present method, using the moving least squares approximation, is already smooth enough. The numerical results in these examples show that high rates of convergence for the Sobolev norms ∥·∥0 and ∥·∥1 are achievable, and that the values of the unknown variable and its derivatives are quite accurate.