[oneD] Only evaluate 'V' for strained flames

Cantera · Aug 13, 2023 · f202340 · f202340
1 parent 5152b22
commit f202340
Show file tree

Hide file tree

Showing 2 changed files with 16 additions and 10 deletions.
diff --git a/src/oneD/Boundary1D.cpp b/src/oneD/Boundary1D.cpp
@@ -183,10 +183,6 @@ void Inlet1D::eval(size_t jg, double* xg, double* rg,
         // the inlet, since this is set within the flow domain from the
         // continuity equation.
 
-        // spreading rate. The flow domain sets this to V(0),
-        // so for finite spreading rate subtract m_V0.
-        rb[c_offset_V] -= m_V0;
-
         if (m_flow->doEnergy(0)) {
             // The third flow residual is for T, where it is set to T(0).  Subtract
             // the local temperature to hold the flow T to the inlet T.
@@ -204,6 +200,10 @@ void Inlet1D::eval(size_t jg, double* xg, double* rg,
             // The flow domain sets this to -rho*u. Add mdot to specify the mass
             // flow rate
             rb[c_offset_L] += m_mdot;
+
+            // spreading rate. The flow domain sets this to V(0),
+            // so for finite spreading rate subtract m_V0.
+            rb[c_offset_V] -= m_V0;
         } else {
             rb[c_offset_U] = m_flow->density(0) * xb[c_offset_U] - m_mdot;
             rb[c_offset_L] = xb[c_offset_L];
@@ -217,7 +217,7 @@ void Inlet1D::eval(size_t jg, double* xg, double* rg,
         }
 
     } else {
-        // right inlet
+        // right inlet (should only be used for counter-flow flames)
         // Array elements corresponding to the last point in the flow domain
         double* rb = rg + loc() - m_flow->nComponents();
         rb[c_offset_V] -= m_V0;

diff --git a/src/oneD/StFlow.cpp b/src/oneD/StFlow.cpp
@@ -513,11 +513,16 @@ void StFlow::evalResidual(double* x, double* rsd, int* diag,
             //    \rho dV/dt + \rho u dV/dz + \rho V^2
             //       = d(\mu dV/dz)/dz - lambda
             //-------------------------------------------------
-            rsd[index(c_offset_V,j)]
-            = (shear(x,j) - lambda(x,j) - rho_u(x,j)*dVdz(x,j)
-               - m_rho[j]*V(x,j)*V(x,j))/m_rho[j]
-              - rdt*(V(x,j) - V_prev(j));
-            diag[index(c_offset_V, j)] = 1;
+            if (m_usesLambda) {
+                rsd[index(c_offset_V,j)] =
+                    (shear(x, j) - lambda(x, j) - rho_u(x, j) * dVdz(x, j)
+                    - m_rho[j] * V(x, j) * V(x, j)) / m_rho[j]
+                    - rdt * (V(x, j) - V_prev(j));
+                diag[index(c_offset_V, j)] = 1;
+            } else {
+                rsd[index(c_offset_V, j)] = V(x, j);
+                diag[index(c_offset_V, j)] = 0;
+            }
 
             //-------------------------------------------------
             //    Species equations
@@ -1010,6 +1015,7 @@ void StFlow::evalRightBoundary(double* x, double* rsd, int* diag, double rdt)
     // and T, and zero diffusive flux for all species.
 
     rsd[index(c_offset_V,j)] = V(x,j);
+    diag[index(c_offset_V,j)] = 0;
     double sum = 0.0;
     // set residual of poisson's equ to zero
     rsd[index(c_offset_E, j)] = x[index(c_offset_E, j)];