Much of the algorithm for relaxing the element potentials
to equilibrium employs the Brinkley algorithm (const T, P)
(p. 129 Smith and Missen). This algorithm does better at
relaxing the system towards steady state from a bad initial
guess.
The old algorithm now does the last little bit and it does
the calculation for fixed pairs other than const T, P, using
the Brinkley algorithm's initial guess.
Modified the damping algorithm in the old method. It now doesn't
create false convergence failures when it is close to the
end of the calculation.
The new solver has been tested against a matrix of initial conditions
for the gri30 mechanism and a 200 species jp8 mechanism.
Yet to do:
This routine should work but has not been checked out against:
phases with electric charge
nonideal phases.
This routine should probably not yet work with
molality based phases.
More work may have to be done with the Brinkley algorithm for
problems which suffer from stoichiometric constraints. A full
Gauss-Jordon elimination step should probably be added to
flag cases where there are effective stoichiometric constraints.
it sums to one. Previously, it was set uniformly to zero. The mass
fraction of species 0 is set to 1.
This was noticed by a user, and may have resulted in seg faults before
the ThermoPhase object was properly initialized.
of progress. This bug was reported in the Cantera newsgroup.
The fix consists of replacing
Rc = exp ( mlogA + blogT - Ea / RT)
with
Rc = A * exp (b logT - Ea/RT)
Therefore, the log of A never has to be taken.
Note, this also allows for a zero rate of progress for a reaction.
Consequences:
The permissible range of the rate constant values may be altered.
Numerical roundoff differences will occur, since the order of
operations in key rate constant evalulations has been changed.
Speed of the calculation may be affected. I have no idea which way.
The fix should be considered as provisional. There are different pathways
for fixing this. This is merely the simplest.