research/cantera
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

[Reactor] Refactor handling of sensitivity parameters in evalEqs

Browse source
This commit is contained in:
Ray Speth 2013-07-30 15:57:48 +00:00
parent 957c88fcff
commit 0295f975e0
5 changed files with 53 additions and 130 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

6
include/cantera/zeroD/Reactor.h
View file

@ -145,6 +145,12 @@ public:
virtual size_t componentIndex(const std::string& nm) const; virtual size_t componentIndex(const std::string& nm) const;
protected: protected:
//! Set reaction rate multipliers based on the sensitivity variables in
//! *params*.
virtual void applySensitivity(double* params);
//! Reset the reaction rate multipliers
virtual void resetSensitivity(double* params);
//! Pointer to the homogeneous Kinetics object that handles the reactions //! Pointer to the homogeneous Kinetics object that handles the reactions
Kinetics* m_kin; Kinetics* m_kin;

35
src/zeroD/ConstPressureReactor.cpp
View file

@ -124,25 +124,7 @@ void ConstPressureReactor::evalEqs(doublereal time, doublereal* y,
m_thermo->restoreState(m_state); m_thermo->restoreState(m_state);
Kinetics* kin; Kinetics* kin;
size_t npar, ploc; applySensitivity(params);
double mult;
// process sensitivity parameters
if (params) {
npar = m_pnum.size();
for (size_t n = 0; n < npar; n++) {
mult = m_kin->multiplier(m_pnum[n]);
m_kin->setMultiplier(m_pnum[n], mult*params[n]);
}
ploc = npar;
for (size_t m = 0; m < m_nwalls; m++) {
if (m_nsens_wall[m] > 0) {
m_wall[m]->setSensitivityParameters(m_lr[m], params + ploc);
ploc += m_nsens_wall[m];
}
}
}
m_Q = 0.0; m_Q = 0.0;
@ -236,20 +218,7 @@ void ConstPressureReactor::evalEqs(doublereal time, doublereal* y,
} }
// reset sensitivity parameters // reset sensitivity parameters
if (params) { resetSensitivity(params);
npar = m_pnum.size();
for (size_t n = 0; n < npar; n++) {
mult = m_kin->multiplier(m_pnum[n]);
m_kin->setMultiplier(m_pnum[n], mult/params[n]);
}
ploc = npar;
for (size_t m = 0; m < m_nwalls; m++) {
if (m_nsens_wall[m] > 0) {
m_wall[m]->resetSensitivityParameters(m_lr[m]);
ploc += m_nsens_wall[m];
}
}
}
} }
size_t ConstPressureReactor::componentIndex(const string& nm) const size_t ConstPressureReactor::componentIndex(const string& nm) const

37
src/zeroD/IdealGasConstPressureReactor.cpp
View file

@ -129,28 +129,9 @@ void IdealGasConstPressureReactor::evalEqs(doublereal time, doublereal* y,
{ {
size_t nk; size_t nk;
m_thermo->restoreState(m_state); m_thermo->restoreState(m_state);
applySensitivity(params);
Kinetics* kin; Kinetics* kin;
size_t npar, ploc;
double mult;
// process sensitivity parameters
if (params) {
npar = m_pnum.size();
for (size_t n = 0; n < npar; n++) {
mult = m_kin->multiplier(m_pnum[n]);
m_kin->setMultiplier(m_pnum[n], mult*params[n]);
}
ploc = npar;
for (size_t m = 0; m < m_nwalls; m++) {
if (m_nsens_wall[m] > 0) {
m_wall[m]->setSensitivityParameters(m_lr[m], params + ploc);
ploc += m_nsens_wall[m];
}
}
}
m_Q = 0.0; m_Q = 0.0;
// compute wall terms // compute wall terms
@ -245,21 +226,7 @@ void IdealGasConstPressureReactor::evalEqs(doublereal time, doublereal* y,
ydot[1] = 0.0; ydot[1] = 0.0;
} }
// reset sensitivity parameters resetSensitivity(params);
if (params) {
npar = m_pnum.size();
for (size_t n = 0; n < npar; n++) {
mult = m_kin->multiplier(m_pnum[n]);
m_kin->setMultiplier(m_pnum[n], mult/params[n]);
}
ploc = npar;
for (size_t m = 0; m < m_nwalls; m++) {
if (m_nsens_wall[m] > 0) {
m_wall[m]->resetSensitivityParameters(m_lr[m]);
ploc += m_nsens_wall[m];
}
}
}
} }
size_t IdealGasConstPressureReactor::componentIndex(const string& nm) const size_t IdealGasConstPressureReactor::componentIndex(const string& nm) const

33
src/zeroD/IdealGasReactor.cpp
View file

@ -138,22 +138,7 @@ void IdealGasReactor::evalEqs(doublereal time, doublereal* y,
doublereal* ydot, doublereal* params) doublereal* ydot, doublereal* params)
{ {
m_thermo->restoreState(m_state); m_thermo->restoreState(m_state);
applySensitivity(params);
// process sensitivity parameters
if (params) {
size_t npar = m_pnum.size();
for (size_t n = 0; n < npar; n++) {
double mult = m_kin->multiplier(m_pnum[n]);
m_kin->setMultiplier(m_pnum[n], mult*params[n]);
}
size_t ploc = npar;
for (size_t m = 0; m < m_nwalls; m++) {
if (m_nsens_wall[m] > 0) {
m_wall[m]->setSensitivityParameters(m_lr[m], params + ploc);
ploc += m_nsens_wall[m];
}
}
}
m_vdot = 0.0; m_vdot = 0.0;
m_Q = 0.0; m_Q = 0.0;
@ -261,21 +246,7 @@ void IdealGasReactor::evalEqs(doublereal time, doublereal* y,
"ydot[" + int2str(i) + "] is not finite"); "ydot[" + int2str(i) + "] is not finite");
} }
// reset sensitivity parameters resetSensitivity(params);
if (params) {
size_t npar = m_pnum.size();
for (size_t n = 0; n < npar; n++) {
double mult = m_kin->multiplier(m_pnum[n]);
m_kin->setMultiplier(m_pnum[n], mult/params[n]);
}
size_t ploc = npar;
for (size_t m = 0; m < m_nwalls; m++) {
if (m_nsens_wall[m] > 0) {
m_wall[m]->resetSensitivityParameters(m_lr[m]);
ploc += m_nsens_wall[m];
}
}
}
} }
size_t IdealGasReactor::componentIndex(const string& nm) const size_t IdealGasReactor::componentIndex(const string& nm) const

72
src/zeroD/Reactor.cpp
View file

@ -180,22 +180,7 @@ void Reactor::evalEqs(doublereal time, doublereal* y,
doublereal* ydot, doublereal* params) doublereal* ydot, doublereal* params)
{ {
m_thermo->restoreState(m_state); m_thermo->restoreState(m_state);
applySensitivity(params);
// process sensitivity parameters
if (params) {
size_t npar = m_pnum.size();
for (size_t n = 0; n < npar; n++) {
double mult = m_kin->multiplier(m_pnum[n]);
m_kin->setMultiplier(m_pnum[n], mult*params[n]);
}
size_t ploc = npar;
for (size_t m = 0; m < m_nwalls; m++) {
if (m_nsens_wall[m] > 0) {
m_wall[m]->setSensitivityParameters(m_lr[m], params + ploc);
ploc += m_nsens_wall[m];
}
}
}
m_vdot = 0.0; m_vdot = 0.0;
m_Q = 0.0; m_Q = 0.0;
@ -306,21 +291,7 @@ void Reactor::evalEqs(doublereal time, doublereal* y,
"ydot[" + int2str(i) + "] is not finite"); "ydot[" + int2str(i) + "] is not finite");
} }
// reset sensitivity parameters resetSensitivity(params);
if (params) {
size_t npar = m_pnum.size();
for (size_t n = 0; n < npar; n++) {
double mult = m_kin->multiplier(m_pnum[n]);
m_kin->setMultiplier(m_pnum[n], mult/params[n]);
}
size_t ploc = npar;
for (size_t m = 0; m < m_nwalls; m++) {
if (m_nsens_wall[m] > 0) {
m_wall[m]->resetSensitivityParameters(m_lr[m]);
ploc += m_nsens_wall[m];
}
}
}
} }
void Reactor::addSensitivityReaction(size_t rxn) void Reactor::addSensitivityReaction(size_t rxn)
@ -383,4 +354,43 @@ size_t Reactor::componentIndex(const string& nm) const
return npos; return npos;
} }
void Reactor::applySensitivity(double* params)
{
if (!params) {
return;
}
size_t npar = m_pnum.size();
for (size_t n = 0; n < npar; n++) {
double mult = m_kin->multiplier(m_pnum[n]);
m_kin->setMultiplier(m_pnum[n], mult*params[n]);
}
size_t ploc = npar;
for (size_t m = 0; m < m_nwalls; m++) {
if (m_nsens_wall[m] > 0) {
m_wall[m]->setSensitivityParameters(m_lr[m], params + ploc);
ploc += m_nsens_wall[m];
}
}
}
void Reactor::resetSensitivity(double* params)
{
if (!params) {
return;
}
size_t npar = m_pnum.size();
for (size_t n = 0; n < npar; n++) {
double mult = m_kin->multiplier(m_pnum[n]);
m_kin->setMultiplier(m_pnum[n], mult/params[n]);
}
size_t ploc = npar;
for (size_t m = 0; m < m_nwalls; m++) {
if (m_nsens_wall[m] > 0) {
m_wall[m]->resetSensitivityParameters(m_lr[m]);
ploc += m_nsens_wall[m];
}
}
}
} }