Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Preconditioned solver improvements #1377

Merged
merged 16 commits into from Aug 30, 2022
Merged

Conversation

speth
Copy link
Member

@speth speth commented Aug 26, 2022

Short version: This restores the performance of the preconditioned solver for reactor networks to what it was like before #1010 was rebased to resolve merge conflicts with the updated Reactor API from #1003.

Changes proposed in this pull request

  • Expand the set of CVODES solver statistics collected and simplify the interface to a single solverStats() function (or solver_stats property in Python.
  • Enable the selection of which side of the system to applying the preconditioner, and make preconditioning on the right the default -- I find that this often significantly reduces the number of time steps needed.
  • Simplify and generalize the conversion of Eigen sparse matrices to Python
  • Make the reactor Jacobian available in Python
  • Fix some issues with the calculation of calculation of derivatives with respect to temperature (done by finite difference in the Jacobian calculation) -- the wrong arguments were being passed into the Reactor::eval(double t, double* lhs, double* rhs) method, treating it like the now-removed Reactor::evalEqs(double t, double* y, double* ydot) method. This resulted in wildly wrong values for these derivatives that caused poor performance for the preconditioned integrator.
  • Fix missing implementation of IdealGasMoleReactor::componentName
  • Make a single-reactor finite difference Jacobian available in both C++ and Python, to aid in verification and debugging.
  • Add missing terms to the energy equation derivatives related to the off-diagonal derivatives of dXj/dni. Add comments to clarify that these terms are deliberately excluded from the species equation to avoid making the Jacobian dense.

If applicable, provide an example illustrating new features this pull request is introducing

Results for a sample of ignition delay calculations before this PR:

vars =   11, steps =  420, time = 0.003, t_ig = 2.04e-03
vars =   54, steps = 5860, time = 0.331, t_ig = 7.93e-01
vars =  121, failed
vars =  202, steps = 1639, time = 0.310, t_ig = 1.54e-02
vars =  231, failed
vars =  494, steps = 2072, time = 0.965, t_ig = 7.34e-02
vars =  655, failed
vars =  875, failed
vars = 1379, steps = 1504, time = 1.447, t_ig = 3.84e-02
vars = 2116, steps = 1702, time = 2.450, t_ig = 2.22e-02
vars = 2650, steps = 9144, time = 29.482, t_ig = 1.74e-02
vars = 3788, failed
vars = 7172, steps = 1878, time = 20.547, t_ig = 5.40e-02

Results for left preconditioning, using this PR:

vars =   11, steps =  420, time = 0.002, t_ig = 2.04e-03
vars =   54, steps = 1106, time = 0.043, t_ig = 7.93e-01
vars =  121, steps = 1444, time = 0.152, t_ig = 5.56e-02
vars =  202, steps = 1100, time = 0.182, t_ig = 1.55e-02
vars =  231, steps = 4156, time = 0.927, t_ig = 3.84e-02
vars =  494, steps = 1088, time = 0.412, t_ig = 7.34e-02
vars =  655, steps = 4190, time = 2.065, t_ig = 2.69e-02
vars =  875, steps = 4926, time = 3.398, t_ig = 1.41e-01
vars = 1379, steps = 1207, time = 0.963, t_ig = 3.84e-02
vars = 2116, steps = 1122, time = 1.488, t_ig = 2.22e-02
vars = 2650, steps = 6480, time = 14.042, t_ig = 1.74e-02
vars = 3788, steps = 6097, time = 16.011, t_ig = 1.54e-02
vars = 7172, steps = 1365, time = 9.583, t_ig = 5.40e-02

Results for right preconditioning, using this PR:

vars =   11, steps =  415, time = 0.002, t_ig = 2.04e-03
vars =   54, steps = 1026, time = 0.044, t_ig = 7.93e-01
vars =  121, steps =  975, time = 0.101, t_ig = 5.56e-02
vars =  202, steps =  918, time = 0.157, t_ig = 1.54e-02
vars =  231, steps = 1150, time = 0.281, t_ig = 3.85e-02
vars =  494, steps =  899, time = 0.380, t_ig = 7.34e-02
vars =  655, steps = 1163, time = 0.635, t_ig = 2.69e-02
vars =  875, steps = 1445, time = 1.055, t_ig = 1.41e-01
vars = 1379, steps =  933, time = 0.808, t_ig = 3.84e-02
vars = 2116, steps =  903, time = 1.212, t_ig = 2.22e-02
vars = 2650, steps = 1296, time = 3.187, t_ig = 1.74e-02
vars = 3788, steps = 1408, time = 4.642, t_ig = 1.54e-02
vars = 7172, steps =  947, time = 7.641, t_ig = 5.40e-02

And for reference, results using the standard dense integrator:

vars =   11, steps =  561, time = 0.002, t_ig = 2.04e-03
vars =   54, steps = 1398, time = 0.030, t_ig = 7.93e-01
vars =  121, steps = 2457, time = 0.268, t_ig = 5.56e-02
vars =  202, steps = 2233, time = 0.513, t_ig = 1.54e-02
vars =  231, steps = 2124, time = 0.532, t_ig = 3.84e-02
vars =  494, steps = 2239, time = 2.105, t_ig = 7.34e-02
vars =  655, steps = 2015, time = 2.680, t_ig = 2.69e-02
vars =  875, steps = 2874, time = 6.401, t_ig = 1.41e-01
vars = 1379, steps = 2212, time = 13.761, t_ig = 3.84e-02
vars = 2116, steps = 2131, time = 31.797, t_ig = 2.22e-02
vars = 2650, steps = 2862, time = 71.558, t_ig = 1.74e-02
vars = 3788, steps = 2510, time = 128.635, t_ig = 1.54e-02
vars = 7172, steps = 2103, time = 531.998, t_ig = 5.40e-02

(and note that this last set is making use of the parallel BLAS/LAPACK implementation available through the macOS Accelerate framework, while the others are all purely single-threaded).

Checklist

  • The pull request includes a clear description of this code change
  • Commit messages have short titles and reference relevant issues
  • Build passes (scons build & scons test) and unit tests address code coverage
  • Style & formatting of contributed code follows contributing guidelines
  • The pull request is ready for review

Using right preconditioning provides a more reliable performance benefit
even after the ignition event, allowing the solver to take reasonable
timesteps and showing a CPU-time benefit even for a single reactor with
the 100-species n-dodecane mechanism.
Finite differences were being calculated incorrectly due to the arguments
to the eval function being misinterpreted.
This is a useful reference for checking analytic or approximate Jacobian
formulations, even if it isn't a particularly good Jacobian to use for
integration generally.
Make exclusion of terms that cause Jacobian to be dense explicit
rather than accidental.

Include full derivatives of species equation in the derivatives of the
energy equation, since this does not affect sparsity.
@codecov
Copy link

codecov bot commented Aug 26, 2022

Codecov Report

Merging #1377 (ac915cf) into main (943ecc6) will increase coverage by 0.05%.
The diff coverage is 79.29%.

@@            Coverage Diff             @@
##             main    #1377      +/-   ##
==========================================
+ Coverage   70.89%   70.95%   +0.05%     
==========================================
  Files         357      357              
  Lines       51500    51532      +32     
  Branches    17198    17252      +54     
==========================================
+ Hits        36510    36562      +52     
+ Misses      12696    12664      -32     
- Partials     2294     2306      +12     
Impacted Files Coverage Δ
include/cantera/cython/kinetics_utils.h 91.80% <ø> (+0.37%) ⬆️
include/cantera/numerics/CVodesIntegrator.h 35.71% <ø> (ø)
...e/cantera/zeroD/IdealGasConstPressureMoleReactor.h 100.00% <ø> (ø)
include/cantera/zeroD/IdealGasMoleReactor.h 100.00% <ø> (ø)
include/cantera/zeroD/Reactor.h 68.18% <0.00%> (+9.09%) ⬆️
include/cantera/zeroD/ReactorNet.h 76.66% <ø> (ø)
include/cantera/numerics/Integrator.h 15.05% <23.07%> (-0.34%) ⬇️
interfaces/cython/cantera/preconditioners.pyx 50.00% <50.00%> (ø)
interfaces/cython/cantera/reactor.pyx 89.23% <50.00%> (+1.16%) ⬆️
src/zeroD/IdealGasMoleReactor.cpp 84.82% <69.04%> (-6.91%) ⬇️
... and 13 more

📣 We’re building smart automated test selection to slash your CI/CD build times. Learn more

@bryanwweber
Copy link
Member

Whoa. Great work @speth!

@ischoegl
Copy link
Member

Nice!

Copy link
Member

@ischoegl ischoegl left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@speth ... I did not look into fine details, as this is something that @anthony-walker may be better prepared to do.

Apart from the minor tweak to make some getters consistent, I was wondering whether the semi-analytical m_kin->netProductionRates_ddT could be leveraged in this context.

include/cantera/numerics/PreconditionerBase.h Outdated Show resolved Hide resolved
interfaces/cython/cantera/preconditioners.pyx Outdated Show resolved Hide resolved
include/cantera/zeroD/Reactor.h Show resolved Hide resolved
test/python/test_reactor.py Show resolved Hide resolved
Copy link
Contributor

@anthony-walker anthony-walker left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@speth Nice work! I think a lot of the renaming and refactoring you have done makes sense.

src/zeroD/Reactor.cpp Show resolved Hide resolved
@@ -9,6 +9,11 @@
namespace Cantera
{

AdaptivePreconditioner::AdaptivePreconditioner()
{
setPreconditionerSide("right");
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It is interesting that "right" preconditioning provides a better result. Does this apply to both constant pressure and constant volume cases?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, right preconditioning seems to generally be better for the constant volume case as well. I have no real insight as to why.

@speth
Copy link
Member Author

speth commented Aug 30, 2022

Apart from the minor tweak to make some getters consistent, I was wondering whether the semi-analytical m_kin->netProductionRates_ddT could be leveraged in this context.

Yes, I think those derivatives could be used, but there are other terms, especially in the energy equation, that would have to be worked out first before you could avoid taking the temperature derivative as a finite difference based on Reactor::eval.

@speth speth requested a review from ischoegl August 30, 2022 13:53
Copy link
Member

@ischoegl ischoegl left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@speth … perfect, thanks!

@speth speth merged commit c7d84a7 into Cantera:main Aug 30, 2022
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

4 participants