Features in the inflaton potential that are traversed in much less than an
e-fold of the expansion can produce observably large non-Gaussianity. In these
models first order corrections to the curvature mode function evolution induce
effects at second order in the slow roll parameters that are generically
greater than ~ 10% and can reach order unity for order unity power spectrum
features. From a complete first order expression in generalized slow-roll, we
devise a computationally efficient method that is as simple to evaluate as the
leading order one and implements consistency relations in a controlled fashion.
This expression matches direct numerical computation for step potential models
of the dominant bispectrum configurations to better than 1% when features are
small and 10% when features are order unity