Wilhoit.fit_to_data might return a bad fit silently

[amarkpayne]

Bug Description

Sometimes, warning messages like the ones below are printed during a call to Wilhoit.fit_to_data

rmg_env/lib/python3.7/site-packages/scipy/optimize/optimize.py:1826: LinAlgWarning: Ill-conditioned matrix (rcond=3.56405e-21): result may not be accurate.
  fu = func(x, *args)

rmgpy/thermo/thermoengine.py:97: LinAlgWarning: Ill-conditioned matrix (rcond=1.08953e-21): result may not be accurate.
  thermo = wilhoit.to_nasa(Tmin=100.0, Tmax=5000.0, Tint=1000.0)

The lines of code in question are here:

RMG-Py/rmgpy/thermo/wilhoit.pyx

Lines 217 to 240 in 447ff5e

	def fit_to_data(self,
	np.ndarray[np.float64_t, ndim=1] Tdata,
	np.ndarray[np.float64_t, ndim=1] Cpdata,
	double Cp0, double CpInf,
	double H298, double S298, double B0=500.0):
	"""
	Fit a Wilhoit model to the data points provided, allowing the
	characteristic temperature `B` to vary so as to improve the fit. This
	procedure requires an optimization, using the ``fminbound`` function
	in the ``scipy.optimize`` module. The data consists of a set
	of heat capacity points `Cpdata` in J/mol*K at a given set of
	temperatures `Tdata` in K, along with the enthalpy `H298` in kJ/mol and
	entropy `S298` in J/mol*K at 298 K. The linearity of the molecule,
	number of vibrational frequencies, and number of internal rotors
	(`linear`, `Nfreq`, and `Nrotors`, respectively) is used to set the
	limits at zero and infinite temperature.
	"""
	self.B = (B0,"K")
	import scipy.optimize
	scipy.optimize.fminbound(self._residual, 300.0, 3000.0, args=(Tdata, Cpdata, Cp0, CpInf, H298, S298))
	return self
	@cython.boundscheck(False)
	@cython.wraparound(False)

As far as I can tell, we don't do any checking here to see if the fit is bad, so unless scipy.optimize.fminbound throws an error (which it doesn't appear to do) this will fail silently. However, fminbound does return a flag that tells you whether the optimization converged or not. It might be worthwhile to use this information (if nothing else print a warning for the specific species, but maybe try re-fitting using another method. I could see why we wouldn't want this to stop an RMG job in its tracks).

How To Reproduce

Run the generate_reactions IPython notebook (once #1735 is merged in or do it on this branch) and you will see the error messages. I will need to do some digging to see which species is causing the issue specifically.

Expected Behavior

If we obtain a bad fit, we should do something about it, and at a minimum, warn the user for the specific species. It is possible that we handle this properly and I am simply missing some code that does this error handling, but it is worth double checking.

[amarkpayne]

The top warning message appears to be related to what @ajocher was seeing in #1726

[kspieks]

I was not able to reproduce the 2 error messages listed at the top of this issue. I switched to ipython branch, made sure RMG-Py and RMG-database were up to date, and then ran generate_reactions.ipynb. The notebook completed, but command F for "LinAlgWarning" or "Ill-conditioned matrix" yielded no results.

Do we still want to try capturing the error flag from scipy.optimize.fminbound (0 if converged, 1 if maximum number of function calls reached) and log it?

[amarkpayne]

Perhaps it depends on the version of the database. I was using the py3_update branch of the database. I'll see if I can replicate it on another computer as well

[sthakris]

Regarding, the issue discussed above. I am running a simulation for large molecule and I get same exact warning and also regarding LinAlgWarning. My model generation is stuck forever, perhaps this issue is causing stiffness in model generation. Can you comment on this? Still no luck with model generation for large molecules like TNT and PETN.

[github-actions]

This issue is being automatically marked as stale because it has not received any interaction in the last 90 days. Please leave a comment if this is still a relevant issue, otherwise it will automatically be closed in 30 days.