Chlorinated solvents are toxic and poorly biodegradable pollutants frequently found in contaminated aquifers. Experimental data of chloroform (CF) aerobic cometabolic biodegradation in a sand column with butane as growth substrate were simulated with a system of non-stationary second-order partial differential equations with non-linear kinetic terms. A MATLAB optimization code based on the Gauss-Newton method and coupled with the Comsol Multiphysics finite elements solver was developed to calibrate the model. For each experimental phase, the best-fit quality was evaluated by an innovative multi-variable model adequacy test. The proposed code solved systems of up to 5 partial differential equations and optimized up to 6 unknown parameters, leading to statistically acceptable best-fits. The optimization of the butane/oxygen pulsed feed led to an 82 % CF biodegradation and to a 0.27 g(CF)/g(butane) transformation yield. When the substrate/pollutant ratio was minimized, the standard model of aerobic cometabolism initially tested required additional terms aimed at taking into account the depletion of reducing energy, in order to attain a statistically acceptable best-fit. This is the first work in which a model of aerobic cometabolism taking into account reducing energy availability was applied to a continuous-flow process. The proposed optimization code can be used for model calibration in a wide range of physical problems described by non-stationary, non-linear partial differential equations, a task that no commercial software can perform. The developed code is made available in the Supplementary Material.