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

[1D] Use c_offset_Y instead of literal value "4"

Browse source
This commit is contained in:
BangShiuh 2016-10-19 17:40:04 -04:00 • committed by Ray Speth
parent 478a62d2af
commit 6d5e28ec46
2 changed files with 23 additions and 22 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

16
src/oneD/StFlow.cpp
View file

@ -14,7 +14,7 @@ namespace Cantera
{
StFlow::StFlow(IdealGasPhase* ph, size_t nsp, size_t points) :
Domain1D(nsp+4, points),
Domain1D(nsp+c_offset_Y, points),
m_press(-1.0),
m_nsp(nsp),
m_thermo(0),
@ -39,7 +39,7 @@ StFlow::StFlow(IdealGasPhase* ph, size_t nsp, size_t points) :
size_t nsp2 = m_thermo->nSpecies();
if (nsp2 != m_nsp) {
m_nsp = nsp2;
Domain1D::resize(m_nsp+4, points);
Domain1D::resize(m_nsp+c_offset_Y, points);
}
// make a local copy of the species molecular weight vector
@ -71,7 +71,7 @@ StFlow::StFlow(IdealGasPhase* ph, size_t nsp, size_t points) :
// mass fraction bounds
for (size_t k = 0; k < m_nsp; k++) {
setBounds(4+k, -1.0e-7, 1.0e5);
setBounds(c_offset_Y+k, -1.0e-7, 1.0e5);
}
//-------------------- grid refinement -------------------------
@ -566,7 +566,7 @@ size_t StFlow::componentIndex(const std::string& name) const
} else if (name=="lambda") {
return 3;
} else {
for (size_t n=4; n<m_nsp+4; n++) {
for (size_t n=c_offset_Y; n<m_nsp+c_offset_Y; n++) {
if (componentName(n)==name) {
return n;
}
@ -670,7 +670,7 @@ void StFlow::restore(const XML_Node& dom, doublereal* soln, int loglevel)
size_t k = m_thermo->speciesIndex(nm);
did_species[k] = 1;
for (size_t j = 0; j < np; j++) {
soln[index(k+4,j)] = x[j];
soln[index(k+c_offset_Y,j)] = x[j];
}
}
} else {
@ -766,7 +766,7 @@ XML_Node& StFlow::save(XML_Node& o, const doublereal* const sol)
addFloatArray(gv,"L",x.size(),x.data(),"N/m^4");
for (size_t k = 0; k < m_nsp; k++) {
soln.getRow(4+k, x.data());
soln.getRow(c_offset_Y+k, x.data());
addFloatArray(gv,m_thermo->speciesName(k),
x.size(),x.data(),"","massFraction");
}
@ -872,7 +872,7 @@ void AxiStagnFlow::evalRightBoundary(doublereal* x, doublereal* rsd,
doublereal sum = 0.0;
for (size_t k = 0; k < m_nsp; k++) {
sum += Y(x,k,j);
rsd[index(k+4,j)] = m_flux(k,j-1) + rho_u(x,j)*Y(x,k,j);
rsd[index(k+c_offset_Y,j)] = m_flux(k,j-1) + rho_u(x,j)*Y(x,k,j);
}
rsd[index(c_offset_Y + rightExcessSpecies(), j)] = 1.0 - sum;
diag[index(c_offset_Y + rightExcessSpecies(), j)] = 0;
@ -925,7 +925,7 @@ void FreeFlame::evalRightBoundary(doublereal* x, doublereal* rsd,
diag[index(c_offset_L, j)] = 0;
for (size_t k = 0; k < m_nsp; k++) {
sum += Y(x,k,j);
rsd[index(k+4,j)] = m_flux(k,j-1) + rho_u(x,j)*Y(x,k,j);
rsd[index(k+c_offset_Y,j)] = m_flux(k,j-1) + rho_u(x,j)*Y(x,k,j);
}
rsd[index(c_offset_Y + rightExcessSpecies(), j)] = 1.0 - sum;
diag[index(c_offset_Y + rightExcessSpecies(), j)] = 0;

29
src/oneD/boundaries1D.cpp
View file

@ -6,6 +6,7 @@
#include "cantera/oneD/Inlet1D.h"
#include "cantera/oneD/OneDim.h"
#include "cantera/base/ctml.h"
#include "cantera/oneD/StFlow.h"
using namespace std;
@ -46,7 +47,7 @@ void Bdry1D::_init(size_t n)
m_left_nv = m_flow_left->nComponents();
m_left_points = m_flow_left->nPoints();
m_left_loc = container().start(m_index-1);
m_left_nsp = m_left_nv - 4;
m_left_nsp = m_left_nv - c_offset_Y;
m_phase_left = &m_flow_left->phase();
} else {
throw CanteraError("Bdry1D::_init",
@ -62,7 +63,7 @@ void Bdry1D::_init(size_t n)
m_flow_right = (StFlow*)&r;
m_right_nv = m_flow_right->nComponents();
m_right_loc = container().start(m_index+1);
m_right_nsp = m_right_nv - 4;
m_right_nsp = m_right_nv - c_offset_Y;
m_phase_right = &m_flow_right->phase();
} else {
throw CanteraError("Bdry1D::_init",
@ -158,7 +159,7 @@ void Inlet1D::init()
}
// components = u, V, T, lambda, + mass fractions
m_nsp = m_flow->nComponents() - 4;
m_nsp = m_flow->nComponents() - c_offset_Y;
m_yin.resize(m_nsp, 0.0);
if (m_xstr != "") {
setMoleFractions(m_xstr);
@ -214,7 +215,7 @@ void Inlet1D::eval(size_t jg, doublereal* xg, doublereal* rg,
// add the convective term to the species residual equations
for (size_t k = 0; k < m_nsp; k++) {
if (k != m_flow_right->leftExcessSpecies()) {
rb[4+k] += x[0]*m_yin[k];
rb[c_offset_Y+k] += x[0]*m_yin[k];
}
}
@ -234,7 +235,7 @@ void Inlet1D::eval(size_t jg, doublereal* xg, doublereal* rg,
rb[0] += x[0]; // u
for (size_t k = 0; k < m_nsp; k++) {
if (k != m_flow_left->rightExcessSpecies()) {
rb[4+k] += x[0]*m_yin[k];
rb[c_offset_Y+k] += x[0]*m_yin[k];
}
}
}
@ -447,7 +448,7 @@ void Outlet1D::eval(size_t jg, doublereal* xg, doublereal* rg, integer* diagg,
double* rb = r + 1;
rb[0] = xb[3];
rb[2] = xb[2] - xb[2 + nc];
for (size_t k = 4; k < nc; k++) {
for (size_t k = c_offset_Y; k < nc; k++) {
rb[k] = xb[k] - xb[k + nc];
}
}
@ -464,8 +465,8 @@ void Outlet1D::eval(size_t jg, doublereal* xg, doublereal* rg, integer* diagg,
}
rb[2] = xb[2] - xb[2 - nc]; // zero T gradient
size_t kSkip = 4 + m_flow_left->rightExcessSpecies();
for (size_t k = 4; k < nc; k++) {
size_t kSkip = c_offset_Y + m_flow_left->rightExcessSpecies();
for (size_t k = c_offset_Y; k < nc; k++) {
if (k != kSkip) {
rb[k] = xb[k] - xb[k - nc]; // zero mass fraction gradient
db[k] = 0;
@ -533,7 +534,7 @@ void OutletRes1D::init()
throw CanteraError("OutletRes1D::init","no flow!");
}
m_nsp = m_flow->nComponents() - 4;
m_nsp = m_flow->nComponents() - c_offset_Y;
m_yres.resize(m_nsp, 0.0);
if (m_xstr != "") {
setMoleFractions(m_xstr);
@ -571,8 +572,8 @@ void OutletRes1D::eval(size_t jg, doublereal* xg, doublereal* rg,
rb[2] = xb[2] - xb[2 + nc];
// specified mass fractions
for (size_t k = 4; k < nc; k++) {
rb[k] = xb[k] - m_yres[k-4];
for (size_t k = c_offset_Y; k < nc; k++) {
rb[k] = xb[k] - m_yres[k-c_offset_Y];
}
}
@ -589,9 +590,9 @@ void OutletRes1D::eval(size_t jg, doublereal* xg, doublereal* rg,
}
rb[2] = xb[2] - m_temp; // zero dT/dz
size_t kSkip = m_flow_left->rightExcessSpecies();
for (size_t k = 4; k < nc; k++) {
for (size_t k = c_offset_Y; k < nc; k++) {
if (k != kSkip) {
rb[k] = xb[k] - m_yres[k-4]; // fixed Y
rb[k] = xb[k] - m_yres[k-c_offset_Y]; // fixed Y
db[k] = 0;
}
}
@ -827,7 +828,7 @@ void ReactingSurf1D::eval(size_t jg, doublereal* xg, doublereal* rg,
size_t nSkip = m_flow_left->rightExcessSpecies();
for (size_t nl = 0; nl < m_left_nsp; nl++) {
if (nl != nSkip) {
rb[4+nl] += m_work[nl]*mwleft[nl];
rb[c_offset_Y+nl] += m_work[nl]*mwleft[nl];
}
}
}