For compact domains with smooth boundaries, we present a surface spline approximation scheme that delivers rates in L-p that are optimal for linear approximation in this setting. This scheme can overcome the boundary effects, observed by Johnson [Constr. Approx., 14 (1998), pp. 429-438], by placing centers with greater density near the boundary. It owes its success to an integral identity employing a minimal number of boundary layer potentials, which, in turn, is derived from the boundary layer potential solution to the Dirichlet problem for the m-fold Laplacian. Furthermore, this integral identity is shown to be the "native space extension" of the target function.