Skip to content

Commit 457acc5

Browse files
jbreue16jbreue16
committed
move parameters to particle model
1 parent f0e9be8 commit 457acc5

17 files changed

+92
-158
lines changed

‎src/libcadet/model/GeneralRateModelDG-LinearSolver.cpp

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -221,7 +221,7 @@ void GeneralRateModelDG::assembleDiscretizedGlobalJacobian(double alpha, Indexer
221221
void GeneralRateModelDG::addTimeDerivativeToJacobianParticleShell(linalg::BandedEigenSparseRowIterator& jac, const Indexer& idxr, double alpha, unsigned int parType)
222222
{
223223
parts::cell::addTimeDerivativeToJacobianParticleShell<linalg::BandedEigenSparseRowIterator, true>(jac, alpha, static_cast<double>(_particle[parType].getPorosity()), _disc.nComp, _disc.nBound + _disc.nComp * parType,
224-
_poreAccessFactor.data() + _disc.nComp * parType, _disc.strideBound[parType], _disc.boundOffset + _disc.nComp * parType, _binding[parType]->reactionQuasiStationarity());
224+
_particle[parType].getPoreAccessfactor(), _disc.strideBound[parType], _disc.boundOffset + _disc.nComp * parType, _binding[parType]->reactionQuasiStationarity());
225225
}
226226

227227
} // namespace model

‎src/libcadet/model/GeneralRateModelDG.cpp

Lines changed: 2 additions & 38 deletions
Original file line numberDiff line numberDiff line change
@@ -375,22 +375,6 @@ bool GeneralRateModelDG::configure(IParameterProvider& paramProvider)
375375
}
376376

377377
// Read vectorial parameters (which may also be section dependent; transport)
378-
_filmDiffusionMode = readAndRegisterMultiplexCompTypeSecParam(paramProvider, _parameters, _filmDiffusion, "FILM_DIFFUSION", _disc.nParType, _disc.nComp, _unitOpIdx);
379-
380-
if ((_filmDiffusion.size() < _disc.nComp * _disc.nParType) || (_filmDiffusion.size() % (_disc.nComp * _disc.nParType) != 0))
381-
throw InvalidParameterException("Number of elements in field FILM_DIFFUSION is not a positive multiple of NCOMP * NPARTYPE (" + std::to_string(_disc.nComp * _disc.nParType) + ")");
382-
383-
if (paramProvider.exists("PORE_ACCESSIBILITY"))
384-
_poreAccessFactorMode = readAndRegisterMultiplexCompTypeSecParam(paramProvider, _parameters, _poreAccessFactor, "PORE_ACCESSIBILITY", _disc.nParType, _disc.nComp, _unitOpIdx);
385-
else
386-
{
387-
_poreAccessFactorMode = MultiplexMode::ComponentType;
388-
_poreAccessFactor = std::vector<cadet::active>(_disc.nComp * _disc.nParType, 1.0);
389-
}
390-
391-
if (_disc.nComp * _disc.nParType != _poreAccessFactor.size())
392-
throw InvalidParameterException("Number of elements in field PORE_ACCESSIBILITY differs from NCOMP * NPARTYPE (" + std::to_string(_disc.nComp * _disc.nParType) + ")");
393-
394378
// Add parameters to map
395379
_parameters[makeParamId(hashString("COL_POROSITY"), _unitOpIdx, CompIndep, ParTypeIndep, BoundStateIndep, ReactionIndep, SectionIndep)] = &_colPorosity;
396380

@@ -552,7 +536,7 @@ void GeneralRateModelDG::notifyDiscontinuousSectionTransition(double t, unsigned
552536

553537
for (int parType = 0; parType < _disc.nParType; parType++)
554538
{
555-
_particle[parType].notifyDiscontinuousSectionTransition(t, secIdx, getSectionDependentSlice(_filmDiffusion, _disc.nComp * _disc.nParType, secIdx), &_poreAccessFactor[0]);
539+
_particle[parType].notifyDiscontinuousSectionTransition(t, secIdx);
556540
}
557541

558542
_disc.curSection = secIdx;
@@ -978,7 +962,7 @@ parts::cell::CellParameters GeneralRateModelDG::makeCellResidualParams(unsigned
978962
_disc.strideBound[parType],
979963
qsReaction,
980964
_particle[parType].getPorosity(),
981-
_poreAccessFactor.data() + _disc.nComp * parType,
965+
_particle[parType].getPoreAccessfactor(),
982966
_binding[parType],
983967
(_dynReaction[parType] && (_dynReaction[parType]->numReactionsCombined() > 0)) ? _dynReaction[parType] : nullptr
984968
};
@@ -1173,10 +1157,6 @@ bool GeneralRateModelDG::setParameter(const ParameterId& pId, double value)
11731157
{
11741158
if (pId.unitOperation == _unitOpIdx)
11751159
{
1176-
if (multiplexCompTypeSecParameterValue(pId, hashString("PORE_ACCESSIBILITY"), _poreAccessFactorMode, _poreAccessFactor, _disc.nParType, _disc.nComp, value, nullptr))
1177-
return true;
1178-
if (multiplexCompTypeSecParameterValue(pId, hashString("FILM_DIFFUSION"), _filmDiffusionMode, _filmDiffusion, _disc.nParType, _disc.nComp, value, nullptr))
1179-
return true;
11801160
const int mpIc = multiplexInitialConditions(pId, value, false);
11811161
if (mpIc > 0)
11821162
return true;
@@ -1268,10 +1248,6 @@ void GeneralRateModelDG::setSensitiveParameterValue(const ParameterId& pId, doub
12681248
{
12691249
if (pId.unitOperation == _unitOpIdx)
12701250
{
1271-
if (multiplexCompTypeSecParameterValue(pId, hashString("PORE_ACCESSIBILITY"), _poreAccessFactorMode, _poreAccessFactor, _disc.nParType, _disc.nComp, value, &_sensParams))
1272-
return;
1273-
if (multiplexCompTypeSecParameterValue(pId, hashString("FILM_DIFFUSION"), _filmDiffusionMode, _filmDiffusion, _disc.nParType, _disc.nComp, value, &_sensParams))
1274-
return;
12751251
if (multiplexInitialConditions(pId, value, true) != 0)
12761252
return;
12771253

@@ -1321,18 +1297,6 @@ bool GeneralRateModelDG::setSensitiveParameter(const ParameterId& pId, unsigned
13211297
{
13221298
if (pId.unitOperation == _unitOpIdx)
13231299
{
1324-
if (multiplexCompTypeSecParameterAD(pId, hashString("PORE_ACCESSIBILITY"), _poreAccessFactorMode, _poreAccessFactor, _disc.nParType, _disc.nComp, adDirection, adValue, _sensParams))
1325-
{
1326-
LOG(Debug) << "Found parameter " << pId << ": Dir " << adDirection << " is set to " << adValue;
1327-
return true;
1328-
}
1329-
1330-
if (multiplexCompTypeSecParameterAD(pId, hashString("FILM_DIFFUSION"), _filmDiffusionMode, _filmDiffusion, _disc.nParType, _disc.nComp, adDirection, adValue, _sensParams))
1331-
{
1332-
LOG(Debug) << "Found parameter " << pId << ": Dir " << adDirection << " is set to " << adValue;
1333-
return true;
1334-
}
1335-
13361300
const int mpIc = multiplexInitialConditions(pId, adDirection, adValue);
13371301
if (mpIc > 0)
13381302
{

‎src/libcadet/model/GeneralRateModelDG.hpp

Lines changed: 0 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -311,12 +311,6 @@ class GeneralRateModelDG : public UnitOperationBase
311311
active _colPorosity; //!< Column porosity (external porosity) \f$ \varepsilon_c \f$
312312
std::vector<active> _parTypeVolFrac; //!< Volume fraction of each particle type
313313

314-
// Vectorial parameters
315-
std::vector<active> _filmDiffusion; //!< Film diffusion coefficient \f$ k_f \f$
316-
MultiplexMode _filmDiffusionMode;
317-
std::vector<active> _poreAccessFactor; //!< Pore accessibility factor \f$ F_{\text{acc}} \f$
318-
MultiplexMode _poreAccessFactorMode;
319-
320314
bool _axiallyConstantParTypeVolFrac; //!< Determines whether particle type volume fraction is homogeneous across axial coordinate
321315
bool _analyticJac; //!< Determines whether AD or analytic Jacobians are used
322316
unsigned int _jacobianAdDirs; //!< Number of AD seed vectors required for Jacobian computation

‎src/libcadet/model/LumpedRateModelWithPoresDG-InitialConditions.cpp

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -375,7 +375,7 @@ namespace cadet
375375
_disc.strideBound[type],
376376
_binding[type]->reactionQuasiStationarity(),
377377
_particle[type].getPorosity(),
378-
_poreAccessFactor.data() + _disc.nComp * type,
378+
_particle[type].getPoreAccessfactor(),
379379
_binding[type],
380380
(_dynReaction[type] && (_dynReaction[type]->numReactionsCombined() > 0)) ? _dynReaction[type] : nullptr
381381
};

‎src/libcadet/model/LumpedRateModelWithPoresDG-LinearSolver.cpp

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -202,7 +202,7 @@ void LumpedRateModelWithPoresDG::addTimeDerivativeToJacobianParticleBlock(linalg
202202
// Add derivative with respect to dc_p / dt to Jacobian
203203
jac[0] += alpha;
204204

205-
const double invBetaP = (1.0 - static_cast<double>(_particle[parType].getPorosity())) / (static_cast<double>(_poreAccessFactor[parType * _disc.nComp + comp]) * static_cast<double>(_particle[parType].getPorosity()));
205+
const double invBetaP = (1.0 - static_cast<double>(_particle[parType].getPorosity())) / (static_cast<double>(_particle[parType].getPoreAccessfactor()[comp]) * static_cast<double>(_particle[parType].getPorosity()));
206206

207207
// Add derivative with respect to dq / dt to Jacobian
208208
const int nBound = static_cast<int>(_disc.nBound[parType * _disc.nComp + comp]);

‎src/libcadet/model/LumpedRateModelWithPoresDG.cpp

Lines changed: 2 additions & 40 deletions
Original file line numberDiff line numberDiff line change
@@ -359,17 +359,6 @@ bool LumpedRateModelWithPoresDG::configure(IParameterProvider& paramProvider)
359359
// Read geometry parameters
360360
_colPorosity = paramProvider.getDouble("COL_POROSITY");
361361

362-
// Read vectorial parameters (which may also be section dependent; transport)
363-
_filmDiffusionMode = readAndRegisterMultiplexCompTypeSecParam(paramProvider, _parameters, _filmDiffusion, "FILM_DIFFUSION", _disc.nParType, _disc.nComp, _unitOpIdx);
364-
365-
if (paramProvider.exists("PORE_ACCESSIBILITY"))
366-
_poreAccessFactorMode = readAndRegisterMultiplexCompTypeSecParam(paramProvider, _parameters, _poreAccessFactor, "PORE_ACCESSIBILITY", _disc.nParType, _disc.nComp, _unitOpIdx);
367-
else
368-
{
369-
_poreAccessFactorMode = MultiplexMode::ComponentType;
370-
_poreAccessFactor = std::vector<cadet::active>(_disc.nComp * _disc.nParType, 1.0);
371-
}
372-
373362
// Check whether PAR_TYPE_VOLFRAC is required or not
374363
if ((_disc.nParType > 1) && !paramProvider.exists("PAR_TYPE_VOLFRAC"))
375364
throw InvalidParameterException("The required parameter \"PAR_TYPE_VOLFRAC\" was not found");
@@ -400,11 +389,6 @@ bool LumpedRateModelWithPoresDG::configure(IParameterProvider& paramProvider)
400389
if (_disc.nParType * _disc.nPoints != _parTypeVolFrac.size())
401390
throw InvalidParameterException("Number of elements in field PAR_TYPE_VOLFRAC does not match number of particle types");
402391

403-
if ((_filmDiffusion.size() < _disc.nComp * _disc.nParType) || (_filmDiffusion.size() % (_disc.nComp * _disc.nParType) != 0))
404-
throw InvalidParameterException("Number of elements in field FILM_DIFFUSION is not a positive multiple of NCOMP * NPARTYPE (" + std::to_string(_disc.nComp * _disc.nParType) + ")");
405-
if (_disc.nComp * _disc.nParType != _poreAccessFactor.size())
406-
throw InvalidParameterException("Number of elements in field PORE_ACCESSIBILITY differs from NCOMP * NPARTYPE (" + std::to_string(_disc.nComp * _disc.nParType) + ")");
407-
408392
// Check that particle volume fractions sum to 1.0
409393
for (unsigned int i = 0; i < _disc.nPoints; ++i)
410394
{
@@ -611,7 +595,7 @@ void LumpedRateModelWithPoresDG::notifyDiscontinuousSectionTransition(double t,
611595
_convDispOp.notifyDiscontinuousSectionTransition(t, secIdx, _jacInlet);
612596

613597
for (int parType = 0; parType < _disc.nParType; parType++)
614-
_particle[parType].notifyDiscontinuousSectionTransition(t, secIdx, &_filmDiffusion[0], &_poreAccessFactor[0]);
598+
_particle[parType].notifyDiscontinuousSectionTransition(t, secIdx);
615599
}
616600

617601
void LumpedRateModelWithPoresDG::setFlowRates(active const* in, active const* out) CADET_NOEXCEPT
@@ -1185,7 +1169,7 @@ void LumpedRateModelWithPoresDG::multiplyWithDerivativeJacobian(const Simulation
11851169
// Add derivative with respect to dc_p / dt to Jacobian
11861170
localRet[comp] = localSdot[comp];
11871171

1188-
const double invBetaP = (1.0 - static_cast<double>(_particle[type].getPorosity())) / (static_cast<double>(_poreAccessFactor[type * _disc.nComp + comp]) * static_cast<double>(_particle[type].getPorosity()));
1172+
const double invBetaP = (1.0 - static_cast<double>(_particle[type].getPorosity())) / (static_cast<double>(_particle[type].getPoreAccessfactor()[comp]) * static_cast<double>(_particle[type].getPorosity()));
11891173

11901174
// Add derivative with respect to dq / dt to Jacobian (normal equations)
11911175
for (unsigned int i = 0; i < nBound[comp]; ++i)
@@ -1275,10 +1259,6 @@ bool LumpedRateModelWithPoresDG::setParameter(const ParameterId& pId, double val
12751259

12761260
return true;
12771261
}
1278-
if (multiplexCompTypeSecParameterValue(pId, hashString("FILM_DIFFUSION"), _filmDiffusionMode, _filmDiffusion, _disc.nParType, _disc.nComp, value, nullptr))
1279-
return true;
1280-
if (multiplexCompTypeSecParameterValue(pId, hashString("PORE_ACCESSIBILITY"), _poreAccessFactorMode, _poreAccessFactor, _disc.nParType, _disc.nComp, value, nullptr))
1281-
return true;
12821262

12831263
const int mpIc = multiplexInitialConditions(pId, value, false);
12841264
if (mpIc > 0)
@@ -1344,12 +1324,6 @@ void LumpedRateModelWithPoresDG::setSensitiveParameterValue(const ParameterId& p
13441324

13451325
return;
13461326
}
1347-
if (multiplexCompTypeSecParameterValue(pId, hashString("FILM_DIFFUSION"), _filmDiffusionMode, _filmDiffusion, _disc.nParType, _disc.nComp, value, &_sensParams))
1348-
return;
1349-
if (multiplexCompTypeSecParameterValue(pId, hashString("PORE_ACCESSIBILITY"), _poreAccessFactorMode, _poreAccessFactor, _disc.nParType, _disc.nComp, value, &_sensParams))
1350-
return;
1351-
if (multiplexInitialConditions(pId, value, true) != 0)
1352-
return;
13531327

13541328
if (_convDispOp.setSensitiveParameterValue(_sensParams, pId, value))
13551329
return;
@@ -1387,18 +1361,6 @@ bool LumpedRateModelWithPoresDG::setSensitiveParameter(const ParameterId& pId, u
13871361
return true;
13881362
}
13891363

1390-
if (multiplexCompTypeSecParameterAD(pId, hashString("FILM_DIFFUSION"), _filmDiffusionMode, _filmDiffusion, _disc.nParType, _disc.nComp, adDirection, adValue, _sensParams))
1391-
{
1392-
LOG(Debug) << "Found parameter " << pId << ": Dir " << adDirection << " is set to " << adValue;
1393-
return true;
1394-
}
1395-
1396-
if (multiplexCompTypeSecParameterAD(pId, hashString("PORE_ACCESSIBILITY"), _poreAccessFactorMode, _poreAccessFactor, _disc.nParType, _disc.nComp, adDirection, adValue, _sensParams))
1397-
{
1398-
LOG(Debug) << "Found parameter " << pId << ": Dir " << adDirection << " is set to " << adValue;
1399-
return true;
1400-
}
1401-
14021364
const int mpIc = multiplexInitialConditions(pId, adDirection, adValue);
14031365
if (mpIc > 0)
14041366
{

‎src/libcadet/model/LumpedRateModelWithPoresDG.hpp

Lines changed: 0 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -280,12 +280,6 @@ class LumpedRateModelWithPoresDG : public UnitOperationBase
280280
active _colPorosity; //!< Column porosity (external porosity) \f$ \varepsilon_c \f$
281281
std::vector<active> _parTypeVolFrac; //!< Volume fraction of each particle type
282282

283-
// Vectorial parameters
284-
std::vector<active> _filmDiffusion; //!< Film diffusion coefficient \f$ k_f \f$
285-
MultiplexMode _filmDiffusionMode;
286-
std::vector<active> _poreAccessFactor; //!< Pore accessibility factor \f$ F_{\text{acc}} \f$
287-
MultiplexMode _poreAccessFactorMode;
288-
289283
bool _axiallyConstantParTypeVolFrac; //!< Determines whether particle type volume fraction is homogeneous across axial coordinate
290284
bool _analyticJac; //!< Determines whether AD or analytic Jacobians are used
291285
unsigned int _jacobianAdDirs; //!< Number of AD seed vectors required for Jacobian computation

‎src/libcadet/model/ParameterMultiplexing.cpp

Lines changed: 0 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -716,10 +716,6 @@ MultiplexMode readAndRegisterSingleTypeMultiplexCompTypeSecParam(IParameterProvi
716716
{
717717
case MultiplexMode::Component:
718718
{
719-
std::vector<active> p(nComp);
720-
std::copy(values.begin(), values.end(), p.begin());
721-
values = std::move(p);
722-
723719
if (parTypeIdx == 0)
724720
{
725721
for (unsigned int s = 0; s < nComp; ++s)

‎src/libcadet/model/particle/GeneralRateParticle.cpp

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -220,9 +220,9 @@ namespace model
220220
return parTransportConfigSuccess && bindingConfSuccess && dynReactionConfSuccess;
221221
}
222222

223-
bool GeneralRateParticle::notifyDiscontinuousSectionTransition(double t, unsigned int secIdx, active const* const filmDiff, active const* const poreAccessFactor)
223+
bool GeneralRateParticle::notifyDiscontinuousSectionTransition(double t, unsigned int secIdx)
224224
{
225-
return _parDiffOp->notifyDiscontinuousSectionTransition(t, secIdx, filmDiff, poreAccessFactor);
225+
return _parDiffOp->notifyDiscontinuousSectionTransition(t, secIdx);
226226
}
227227

228228
int GeneralRateParticle::writeParticleCoordinates(double* coords) const

‎src/libcadet/model/particle/GeneralRateParticle.hpp

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -92,7 +92,7 @@ namespace parts
9292

9393
parts::cell::CellParameters makeCellResidualParams(int const* qsReaction, unsigned int const* nBound) const override;
9494

95-
bool notifyDiscontinuousSectionTransition(double t, unsigned int secIdx, active const* const filmDiff, active const* const poreAccessFactor) override;
95+
bool notifyDiscontinuousSectionTransition(double t, unsigned int secIdx) override;
9696

9797
int residual(double t, unsigned int secIdx, double const* yPar, double const* yBulk, double const* yDotPar, double* resPar, double* resBulk, columnPackingParameters packing, linalg::BandedEigenSparseRowIterator& jacIt, LinearBufferAllocator tlmAlloc, WithoutParamSensitivity) override;
9898
int residual(double t, unsigned int secIdx, double const* yPar, double const* yBulk, double const* yDotPar, active* resPar, active* resBulk, columnPackingParameters packing, linalg::BandedEigenSparseRowIterator& jacIt, LinearBufferAllocator tlmAlloc, WithParamSensitivity) override;
@@ -110,6 +110,7 @@ namespace parts
110110

111111
inline const active& getPorosity() const CADET_NOEXCEPT override { return _parDiffOp->getPorosity(); }
112112
inline const active* getPoreAccessfactor() const CADET_NOEXCEPT override { return _parDiffOp->getPoreAccessfactor(); }
113+
inline const active* getFilmDiffusion(const unsigned int secIdx) const CADET_NOEXCEPT { return _parDiffOp->getFilmDiffusion(secIdx); }
113114

114115
inline int nDiscPoints() const CADET_NOEXCEPT override { return _parDiffOp->nDiscPoints(); }
115116
inline int strideParBlock() const CADET_NOEXCEPT override { return nDiscPoints() * stridePoint(); }

0 commit comments

Comments
 (0)