  In this paper, we propose a class of nonlocal models to approximate the
Poisson model on manifolds with homogeneous Neumann boundary condition, where
the manifolds are assumed to be embedded in high dimensional Euclid spaces. In
comparison to the existing nonlocal approximation of Poisson models with
Neumann boundary, we optimize the truncation error of model by adding an
augmented term along the $2\delta$ layer of boundary, with $2\delta$ be the
nonlocal interaction horizon. Such term is formulated by the integration of the
second order normal derivative of solution through the boundary, while the
second order normal derivative is expressed as the difference between the
interior Laplacian and the boundary Laplacian. The concentration of our paper
is on the construction of nonlocal model, the well-posedness of model, and its
second-order convergence rate to its local counterpart. The localization rate
of our nonlocal model is currently optimal among all related works even for the
case of high dimensional Euclid spaces.
