diff --git a/include/cantera/oneD/Sim1D.h b/include/cantera/oneD/Sim1D.h index 083f189c2..a28e06a23 100644 --- a/include/cantera/oneD/Sim1D.h +++ b/include/cantera/oneD/Sim1D.h @@ -101,8 +101,16 @@ public: void setRefineCriteria(int dom = -1, doublereal ratio = 10.0, doublereal slope = 0.8, doublereal curve = 0.8, doublereal prune = -0.1); void setMaxGridPoints(int dom = -1, int npoints = 300); - + void restore(std::string fname, std::string id); + //! Set the minimum grid spacing in the specified domain(s). + /*! + * @param dom Domain index. If dom == -1, the specified spacing + is applied to all domains. + @param gridmin The minimum allowable grid spacing [m] + */ + void setGridMin(int dom, double gridmin); + void getInitialSoln(); void setSolution(const doublereal* soln) { diff --git a/include/cantera/oneD/refine.h b/include/cantera/oneD/refine.h index 2cc47adeb..3491ef9b8 100644 --- a/include/cantera/oneD/refine.h +++ b/include/cantera/oneD/refine.h @@ -29,6 +29,18 @@ public: void setMaxPoints(int npmax) { m_npmax = npmax; } + + //! Set the minimum allowable spacing between adjacent grid points [m]. + void setGridMin(double gridmin) { + m_gridmin = gridmin; + } + + //! Returns the the minimum allowable spacing between adjacent + //! grid points [m]. + double gridMin() const { + return m_gridmin; + } + int analyze(size_t n, const doublereal* z, const doublereal* x); int getNewGrid(int n, const doublereal* z, int nn, doublereal* znew); //int getNewSoln(int n, const doublereal* x, doublereal* xnew); @@ -67,6 +79,7 @@ protected: Domain1D* m_domain; size_t m_nv, m_npmax; doublereal m_thresh; + doublereal m_gridmin; //!< minimum grid spacing [m] }; diff --git a/src/oneD/Sim1D.cpp b/src/oneD/Sim1D.cpp index 61904f31d..b37b620e0 100644 --- a/src/oneD/Sim1D.cpp +++ b/src/oneD/Sim1D.cpp @@ -578,6 +578,20 @@ void Sim1D::setRefineCriteria(int dom, doublereal ratio, } } +void Sim1D::setGridMin(int dom, double gridmin) +{ + if (dom >= 0) { + Refiner& r = domain(dom).refiner(); + r.setGridMin(gridmin); + } else { + for (size_t n = 0; n < m_nd; n++) { + Refiner& r = domain(n).refiner(); + r.setGridMin(gridmin); + } + } +} + + void Sim1D::setMaxGridPoints(int dom, int npoints) { if (dom >= 0) { diff --git a/src/oneD/refine.cpp b/src/oneD/refine.cpp index 765480ff0..ded8fdbf0 100644 --- a/src/oneD/refine.cpp +++ b/src/oneD/refine.cpp @@ -40,7 +40,8 @@ static doublereal eps() Refiner::Refiner(Domain1D& domain) : m_ratio(10.0), m_slope(0.8), m_curve(0.8), m_prune(-0.001), - m_min_range(0.01), m_domain(&domain), m_npmax(3000) + m_min_range(0.01), m_domain(&domain), m_npmax(3000), + m_gridmin(5e-6) { m_nv = m_domain->nComponents(); m_active.resize(m_nv, true); @@ -129,7 +130,7 @@ int Refiner::analyze(size_t n, const doublereal* z, dmax = m_slope*(vmax - vmin) + m_thresh; for (j = 0; j < n-1; j++) { r = fabs(v[j+1] - v[j])/dmax; - if (r > 1.0) { + if (r > 1.0 && dz[j] >= 2 * m_gridmin) { m_loc[j] = 1; m_c[name] = 1; //if (int(m_loc.size()) + n > m_npmax) goto done; @@ -161,7 +162,8 @@ int Refiner::analyze(size_t n, const doublereal* z, dmax = m_curve*(smax - smin); // + 0.5*m_curve*(smax + smin); for (j = 0; j < n-2; j++) { r = fabs(s[j+1] - s[j]) / (dmax + m_thresh/dz[j]); - if (r > 1.0) { + if (r > 1.0 && dz[j] >= 2 * m_gridmin && + dz[j+1] >= 2 * m_gridmin) { m_c[name] = 1; m_loc[j] = 1; m_loc[j+1] = 1;