diff --git a/src/oneD/StFlow.cpp b/src/oneD/StFlow.cpp index aafaed5a0..02fb636ca 100644 --- a/src/oneD/StFlow.cpp +++ b/src/oneD/StFlow.cpp @@ -1054,24 +1054,23 @@ void StFlow::evalContinuity(size_t j, double* x, double* rsd, int* diag, double -(rho_u(x,j+1) - rho_u(x,j))/m_dz[j] -(density(j+1)*V(x,j+1) + density(j)*V(x,j)); } else if (m_isFree) { + // terms involving V are zero as V=0 by definition if (grid(j) > m_zfixed) { rsd[index(c_offset_U,j)] = - - (rho_u(x,j) - rho_u(x,j-1))/m_dz[j-1] - - (density(j-1)*V(x,j-1) + density(j)*V(x,j)); + - (rho_u(x,j) - rho_u(x,j-1))/m_dz[j-1]; } else if (grid(j) == m_zfixed) { if (m_do_energy[j]) { rsd[index(c_offset_U,j)] = (T(x,j) - m_tfixed); } else { rsd[index(c_offset_U,j)] = (rho_u(x,j) - - m_rho[0]*0.3); + - m_rho[0]*0.3); // why 0.3? } } else if (grid(j) < m_zfixed) { rsd[index(c_offset_U,j)] = - - (rho_u(x,j+1) - rho_u(x,j))/m_dz[j] - - (density(j+1)*V(x,j+1) + density(j)*V(x,j)); + - (rho_u(x,j+1) - rho_u(x,j))/m_dz[j]; } } else { - // unstrained + // unstrained with fixed mass flow rate rsd[index(c_offset_U, j)] = rho_u(x, j) - rho_u(x, j - 1); } }