Extensible jacobian and network jacobian interfaces#1634
Extensible jacobian and network jacobian interfaces#1634anthony-walker wants to merge 13 commits intoCantera:mainfrom
Conversation
anthony-walker
force-pushed
the
ext_prec
branch
from
eb36c56 to
205567a
Compare
Codecov Report❌ Patch coverage is Additional details and impacted files@@ Coverage Diff @@
## main #1634 +/- ##
==========================================
+ Coverage 76.75% 76.78% +0.03%
==========================================
Files 457 458 +1
Lines 53744 53929 +185
Branches 9122 9151 +29
==========================================
+ Hits 41250 41411 +161
- Misses 9386 9400 +14
- Partials 3108 3118 +10 ☔ View full report in Codecov by Sentry. 🚀 New features to boost your workflow:
|
anthony-walker
force-pushed
the
ext_prec
branch
3 times, most recently
from
1c005f9 to
ab83aa3
Compare
ischoegl
left a comment
ischoegl
left a comment
There was a problem hiding this comment.
While I don't have enough experience in extensible things for a complete review, I noticed that a lot of test are failing with:
*******************************************************************************
NotImplementedError thrown by MassFlowController::buildReactorJacobian:
Not implemented.
*******************************************************************************
which doesn't appear to be .NET related.
anthony-walker
force-pushed
the
ext_prec
branch
2 times, most recently
from
982eb65 to
96bd8af
Compare
|
Before you get any further, I just wanted to suggest undoing the merge commit and rebasing onto the current |
anthony-walker
force-pushed
the
ext_prec
branch
from
a78547c to
6f61796
Compare
anthony-walker
force-pushed
the
ext_prec
branch
2 times, most recently
from
4857e53 to
84ff012
Compare
anthony-walker
force-pushed
the
ext_prec
branch
from
decbbac to
feb862c
Compare
anthony-walker
force-pushed
the
ext_prec
branch
from
feb862c to
a525724
Compare
anthony-walker
force-pushed
the
ext_prec
branch
from
a525724 to
9a2d1b2
Compare
speth
left a comment
speth
left a comment
There was a problem hiding this comment.
Thanks for putting this together, @anthony-walker. I think it helps in laying out the necessary pieces for handling the cross-reactor terms in the Jacobian -- in particular, the mechanics surrounding Wall::buildNetworkJacobian seem like a good basis to keep building on.
I think there is some room for simplification of this structure, which I've left specific comments on. However, I think there is a significant issue with sorting out exactly which derivatives need to be computed and what the right formulas are for the different cases. The challenge I see us needing to resolve to generalize this across the different reactor types is that what derivative is needed for a single interaction (for example, wall heat transfer) depends on both reactor types and what variables they use to represent each term, and we need a way of selecting the correct formula for each Jacobian element based on that pair of reactor types.
| // derivative of temperature transformed by ideal gas law | ||
| vector<double> moles(m_nsp); | ||
| getMoles(moles.data()); | ||
| double dTdni = pressure() * m_vol / GasConstant / std::accumulate(moles.begin(), moles.end(), 0.0); |
There was a problem hiding this comment.
Like in IdealGasMoleReactor,
test/python/test_reactor.py
Outdated
| # check for values | ||
| for i in range(gas1.n_species): | ||
| assert np.isclose(jac[0, i + 1 + r1.n_vars], - U * A * gas2.T) | ||
| assert np.isclose(jac[r1.n_vars, i + 1], U * A * gas1.T) | ||
| if (i + 1 != 2): | ||
| assert np.isclose(jac[0, i + 1], U * A * gas1.T, rtol=1e-2) | ||
| assert np.isclose(jac[r1.n_vars, i + 1 + r1.n_vars], - U * A * gas2.T, rtol=1e-2) |
There was a problem hiding this comment.
I think you should check against the finite difference Jacobian. This test is just checking against the same incorrect formulas in Python as what's in the C++ implementation.
test/python/test_reactor.py
Outdated
| gas2 = ct.Solution("h2o2.yaml", "ohmech") | ||
| gas2.TPX = 300, ct.one_atm, "O2:1.0" | ||
| r2 = self.reactorClass(gas2) | ||
| r2.chemistry_enabled = False |
There was a problem hiding this comment.
The analytical Jacobian doesn't currently respect the chemistry_enabled flag, so there are a lot of non-zero elements here that should be zero based on the setting of this flag.
|
For an example of where the current derivatives are clearly not correct, consider the following example, adapted from one of your new tests: import cantera as ct
import numpy as np
# create first reactor
gas1 = ct.Solution("h2o2.yaml", "ohmech")
gas1.TPX = 600, ct.one_atm, "O2:1.0"
r1 = ct.IdealGasMoleReactor(gas1)
r1.name = 'R1'
gas1.set_multiplier(0.0)
# create second reactor
gas2 = ct.Solution("h2o2.yaml", "ohmech")
gas2.TPX = 300, ct.one_atm, "O2:1.0"
gas2.set_multiplier(0.0)
r2 = ct.IdealGasMoleReactor(gas2)
r2.name = 'R2'
# create wall
U = 500.0
A = 3.0
w = ct.Wall(r1, r2, U=U, A=A)
net = ct.ReactorNet([r1, r2])
def print_jac(J):
for i,j in zip(*np.nonzero(J)):
name_i = net.component_name(i).replace('temperature', 'T')
name_j = net.component_name(j).replace('temperature', 'T')
net.component_name(i).replace('temperature', 'T')
#print(f"J[{i: 3d}, {j: 3d}] = {J[i,j]:.3e}")
print(f"J[{name_i:8s}, {name_j:8s}] = {J[i,j]:.3e}")
print('Semi-analytical Jacobian:')
print_jac(net.jacobian)
print('\nFinite difference Jacobian:')
print_jac(net.finite_difference_jacobian)which outputs: |
|
@speth apologies for the delay, I have made some updates but still need to do some of the larger changes. It is still on my radar. |
|
@speth pushed up some more changes, I will have to put a bit of thought into some of the other changes. |
anthony-walker
force-pushed
the
ext_prec
branch
from
cf394c3 to
3b892b6
Compare
speth
force-pushed
the
ext_prec
branch
2 times, most recently
from
5f3a0c1 to
0a9d2b4
Compare
anthony-walker
force-pushed
the
ext_prec
branch
from
f3610f9 to
78e660e
Compare
anthony-walker
force-pushed
the
ext_prec
branch
from
78e660e to
42569ef
Compare
It can sometimes be useful to add on to the default evaluation in such a way that the standard before or after delegations cannot handle. A default evaluation function allows for this.
in regards to equations, handling volumes, etc.
3 participants
Changes proposed in this pull request
In this pull request I have setup basic structure for adding contributions to the jacobian from flow devices and walls.
I have also added an extensible jacobian interface to the
ExtensibleReactorsystem and modified the existing jacobian system to pass a vector of triplets as is needed in lieu of creating a large number of separate sparse matrices.@speth @ischoegl some additional tests are needed for the wall and flow contributions but I wanted to get the content officially into a PR.
Checklist
scons build&scons test) and unit tests address code coverage