Predicting the residual saturation of a trapped non-aqueous phase liquid contaminant is critical to estimating the region of contamination, the design of remediation strategies, and risk assessment. Models were developed to predict residual NAPL saturation utilizing optimization and non-linear error functions, consequently allowing for a broader mathematical approach to model development. The input parameters evaluated represent soil and fluid properties: the uniformity coefficient (Cu), the coefficient of gradation (Cc), the capillary number (Nc), the bond number (Nb) and the total trapping number (Nt). Overall, the model that performed best was based on a second-order equation with the independent variables Cu and Nt1 using the sum of the squares of the errors. The nonlinear error function based on a derivative of Marquardt’s Percent Standard Deviation performed best for three other cases.