search

bjodah / chempy

⚗ A package useful for chemistry written in Python
BSD 2-Clause "Simplified" License
536 stars 79 forks source link

Equilibrium root solver fails when some concentrations are 0 #190

Open yakdahl opened 3 years ago

yakdahl commented 3 years ago

Hi! First of all, thanks so much for this package, I really appreciate it!

This is a minor issue: I'm using the Eqsys.root function to solve for equilibrium concentrations in order to analyze a titration curve measured by UV-VIS. I'm trying to estimate extinction coefficients of unknown complexes by simulating these concentrations and noticed that the solver fails when one of my concentrations is 0 i.e. the beginning point of my titration before any addition. While I can work around it, maybe it would be possible to catch this exception before passing the values to the solver and just return the input concentration? Might make more sense than returning a warning and some non-physical concentrations.

Here's an example:

Example_root_fail.txt

bjodah commented 3 years ago

Hi @yakdahl, thanks for letting us know. This is a nice suggestion for improvement.

When a component is completely missing (Cl in this case), we should handle that case by adding code detecting that a component is excluded and solve a smaller problem (with all reactions involving any absent substances removed). In your case, there will be no reactions left, and only the Zn(II) concentrations would be returned.

In general, we cannot guarantee that the solver will always be able to solve an arbitrary problem (due to shortcomings in the numerical solver). So at times there will probably be errors.