@@ -1329,31 +1329,45 @@ int CSTRModel::residualImpl(double t, unsigned int secIdx, StateType const* cons
1329
1329
{
1330
1330
resCMoities[comp] = _MconvMoityBulk.row (MoityIdx) * resC[comp];
1331
1331
MoityIdx++;
1332
+
1333
+ if (wantJac)
1334
+ _jac.native (comp, _nComp + _totalBound) = 0 ; // dF_{Vliquid}/dci = 0
1335
+
1332
1336
}
1333
1337
else if (comp < _nComp - nmoities)
1334
1338
{
1335
1339
resCMoities[comp] += v * flux[comp];
1340
+
1341
+ if (wantJac)
1342
+ _jac.native (comp, _nComp + _totalBound) = static_cast <double >(flux[comp]); // dF/dci = v_liquid
1343
+
1336
1344
}
1345
+
1337
1346
else if (comp > _nComp - nmoities)
1338
1347
{
1339
1348
resCMoities[comp] = v * flux[comp];
1349
+
1350
+ if (wantJac)
1351
+ _dynReactionBulk->analyticJacobianLiquidAdd (t, secIdx, colPos, reinterpret_cast <double const *>(c), -static_cast <double >(v), _jac.row (0 ), subAlloc);
1340
1352
}
1341
1353
}
1354
+
1342
1355
Eigen::Map<Vector<ResidualType, Dynamic>> mapResC (resC, _nComp);
1343
1356
mapResC = resCMoities;
1357
+
1344
1358
}
1345
1359
else
1346
1360
{
1347
1361
for (unsigned int comp = 0 ; comp < _nComp; ++comp)
1348
1362
resC[comp] += v * flux[comp];
1349
- }
1350
-
1351
- if (wantJac)
1352
- {
1353
- for (unsigned int comp = 0 ; comp < _nComp; ++comp)
1354
- _jac.native (comp, _nComp + _totalBound) += static_cast <double >(flux[comp]); // dF/dvliquid = flux
1363
+
1364
+ if (wantJac)
1365
+ {
1366
+ for (unsigned int comp = 0 ; comp < _nComp; ++comp)
1367
+ _jac.native (comp, _nComp + _totalBound) += static_cast <double >(flux[comp]); // dF/dvliquid = flux
1355
1368
1356
- _dynReactionBulk->analyticJacobianLiquidAdd (t, secIdx, colPos, reinterpret_cast <double const *>(c), -static_cast <double >(v), _jac.row (0 ), subAlloc);
1369
+ _dynReactionBulk->analyticJacobianLiquidAdd (t, secIdx, colPos, reinterpret_cast <double const *>(c), -static_cast <double >(v), _jac.row (0 ), subAlloc);
1370
+ }
1357
1371
}
1358
1372
}
1359
1373
0 commit comments