This work explores the numerical translation of the weak or integral solution of nonlinear partial differential equations into a numerically efficient, time-evolving scheme. Specifically, we focus on partial differential equations separable into a quasilinear term and a nonlinear one, with the former defining the Green function of the problem. Utilizing the Green function under a short-time approximation, it becomes possible to derive the integral solution of the problem by breaking it into three integral terms: the propagation of initial conditions and the contributions of the nonlinear and boundary terms. Accordingly, we follow this division to describe and separately analyze the resulting algorithm. To ensure low interpolation error and accurate numerical Green functions, we adapt a piecewise interpolation collocation method to the integral scheme, optimizing the positioning of grid points near the boundary region. At the same time, we employ a second-order quadrature method in time to efficiently implement the nonlinear terms. Validation of both adapted methodologies is conducted by applying them to problems with known analytical solution, as well as to more challenging, norm-preserving problems such as the Burgers equation and the soliton solution of the nonlinear Schrödinger equation. Finally, the boundary term is derived and validated using a series of test cases that cover the range of possible scenarios for boundary problems within the introduced methodology.