From 707c569a3d393de2e0390ca74adf269c6fc85044 Mon Sep 17 00:00:00 2001 From: Ray Speth Date: Mon, 17 Apr 2023 22:45:25 -0400 Subject: [PATCH] [Kinetics] Deprecate reacting phase not in first position Standardize on the behavior implemented by the newSolution and newInterface constructors. --- include/cantera/clib/ct.h | 1 + include/cantera/kinetics/InterfaceKinetics.h | 3 ++- include/cantera/kinetics/Kinetics.h | 3 +++ include/cantera/kinetics/KineticsFactory.h | 7 +++++- interfaces/cython/cantera/kinetics.pyx | 9 +++++++- src/fortran/cantera_funcs.f90 | 4 ++-- src/kinetics/Kinetics.cpp | 5 ++++ src/kinetics/KineticsFactory.cpp | 8 ++++++- test/clib/test_clib.cpp | 2 +- test/clib/test_ctonedim.cpp | 4 ++-- test/clib/test_ctreactor.cpp | 4 ++-- test/general/test_serialization.cpp | 4 ++-- test/kinetics/kineticsFromScratch.cpp | 6 ++--- test/kinetics/rates.cpp | 2 +- test_problems/diamondSurf/runDiamond.cpp | 3 +-- .../surfSolverTest/results2_blessed.txt | 8 +++---- .../surfSolverTest/surfaceSolver2.cpp | 23 ++++++++++--------- 17 files changed, 62 insertions(+), 34 deletions(-) diff --git a/include/cantera/clib/ct.h b/include/cantera/clib/ct.h index b4cb152ab..066ce03ca 100644 --- a/include/cantera/clib/ct.h +++ b/include/cantera/clib/ct.h @@ -122,6 +122,7 @@ extern "C" { CANTERA_CAPI int thermo_setState_Psat(int n, double p, double x); CANTERA_CAPI int thermo_setState_Tsat(int n, double t, double x); + //! @since Starting in Cantera 3.0, the "phasename" argument should be blank CANTERA_CAPI int kin_newFromFile(const char* filename, const char* phasename, int reactingPhase, int neighbor1, int neighbor2, int neighbor3, int neighbor4); diff --git a/include/cantera/kinetics/InterfaceKinetics.h b/include/cantera/kinetics/InterfaceKinetics.h index 8f27cdbb1..50f1fd31e 100644 --- a/include/cantera/kinetics/InterfaceKinetics.h +++ b/include/cantera/kinetics/InterfaceKinetics.h @@ -133,7 +133,8 @@ public: //! Add a thermo phase to the kinetics manager object. /*! * This must be done before the function init() is called or - * before any reactions are input. + * before any reactions are input. The lowest dimensional phase, where reactions + * occur, must be added first. * * This function calls Kinetics::addPhase(). It also sets the following * fields: diff --git a/include/cantera/kinetics/Kinetics.h b/include/cantera/kinetics/Kinetics.h index 20dee5a47..8432f688b 100644 --- a/include/cantera/kinetics/Kinetics.h +++ b/include/cantera/kinetics/Kinetics.h @@ -216,6 +216,9 @@ public: * phase with the smallest spatial dimension (1, 2, or 3) among the list * of phases. If there is more than one, the index of the first one is * returned. For homogeneous mechanisms, the value 0 is returned. + * + * @deprecated Starting in Cantera 3.0, the reacting phase will always be the + * first phase in the InterfaceKinetics object. To be removed after Cantera 3.1. */ size_t reactionPhaseIndex() const { return m_rxnphase; diff --git a/include/cantera/kinetics/KineticsFactory.h b/include/cantera/kinetics/KineticsFactory.h index 703b6d8ea..b55811520 100644 --- a/include/cantera/kinetics/KineticsFactory.h +++ b/include/cantera/kinetics/KineticsFactory.h @@ -82,10 +82,15 @@ unique_ptr newKinetics(const std::vector& phases, * the reactions occur. Searches the Cantera data for this file. * @param phase_name The name of the reacting phase in the input file (that is, the * name of the first phase in the `phases` vector) + * @deprecated The 'phase_name' argument is deprecated and will be removed after + * Cantera 3.0. + * @since Starting with Cantera 3.0, if the reacting phase is not the first item in the + * `phases` vector, a deprecation warning will be issued. In Cantera 3.1, this + * warning will become an error. */ shared_ptr newKinetics(const vector>& phases, const string& filename, - const string& phase_name); + const string& phase_name=""); //! @copydoc newKinetics(const vector>&, const string&, const string&) //! @deprecated To be removed after Cantera 3.0; diff --git a/interfaces/cython/cantera/kinetics.pyx b/interfaces/cython/cantera/kinetics.pyx index fb1afe1c2..bb0f49b47 100644 --- a/interfaces/cython/cantera/kinetics.pyx +++ b/interfaces/cython/cantera/kinetics.pyx @@ -111,7 +111,14 @@ cdef class Kinetics(_SolutionBase): return self.kinetics.nPhases() property reaction_phase_index: - """The index of the phase where the reactions occur.""" + """ + The index of the phase where the reactions occur. + + .. deprecated:: 3.0 + + After Cantera 3.0, the reacting phase is always the first phase associated + with the Kinetics object. This method will be removed after Cantera 3.1. + """ def __get__(self): return self.kinetics.reactionPhaseIndex() diff --git a/src/fortran/cantera_funcs.f90 b/src/fortran/cantera_funcs.f90 index 9496e3fad..5b0d1b897 100644 --- a/src/fortran/cantera_funcs.f90 +++ b/src/fortran/cantera_funcs.f90 @@ -24,7 +24,7 @@ module cantera_funcs self = ctthermo_newFromFile(src) end if - call ctkin_newFromFile(self, src, id) + call ctkin_newFromFile(self, src, '') if (present(loglevel)) then self%tran_id = trans_newdefault(self%thermo_id, loglevel) @@ -52,7 +52,7 @@ module cantera_funcs self%bulk = bulk end if self%gas = gas - call ctkin_newFromFile(self%surf, src, id, gas, bulk) + call ctkin_newFromFile(self%surf, src, '', gas, bulk) ctfunc_importInterface = self return diff --git a/src/kinetics/Kinetics.cpp b/src/kinetics/Kinetics.cpp index 9d4e9ddd9..abe2aea20 100644 --- a/src/kinetics/Kinetics.cpp +++ b/src/kinetics/Kinetics.cpp @@ -596,6 +596,11 @@ void Kinetics::addThermo(shared_ptr thermo) // the phase with lowest dimensionality is assumed to be the // phase/interface at which reactions take place if (thermo->nDim() <= m_mindim) { + if (!m_thermo.empty()) { + warn_deprecated("Kinetics::addThermo", "The reacting (lowest dimensional) " + "phase should be added first. This warning will become an error after " + "Cantera 3.0."); + } m_mindim = thermo->nDim(); m_rxnphase = nPhases(); } diff --git a/src/kinetics/KineticsFactory.cpp b/src/kinetics/KineticsFactory.cpp index 99aab763c..6d643afa7 100644 --- a/src/kinetics/KineticsFactory.cpp +++ b/src/kinetics/KineticsFactory.cpp @@ -137,8 +137,14 @@ shared_ptr newKinetics(const vector>& phases, const string& filename, const string& phase_name) { + if (phase_name != "") { + warn_deprecated("newKinetics", "After Cantera 3.0, the 'phase_name' argument " + "will be removed and an exception will be thrown if the reacting phase is " + "not the first phase in the 'phases' vector."); + } + string reaction_phase = (phase_name.empty()) ? phases.at(0)->name() : phase_name; AnyMap root = AnyMap::fromYamlFile(filename); - AnyMap& phaseNode = root["phases"].getMapWhere("name", phase_name); + AnyMap& phaseNode = root["phases"].getMapWhere("name", reaction_phase); return newKinetics(phases, phaseNode, root); } diff --git a/test/clib/test_clib.cpp b/test/clib/test_clib.cpp index ea7a4ceee..a8382db1c 100644 --- a/test/clib/test_clib.cpp +++ b/test/clib/test_clib.cpp @@ -175,7 +175,7 @@ TEST(ct, thermo) TEST(ct, kinetics) { int thermo = thermo_newFromFile("gri30.yaml", "gri30"); - int kin = kin_newFromFile("gri30.yaml", "gri30", thermo, -1, -1, -1, -1); + int kin = kin_newFromFile("gri30.yaml", "", thermo, -1, -1, -1, -1); ASSERT_GE(kin, 0); size_t nr = kin_nReactions(kin); diff --git a/test/clib/test_ctonedim.cpp b/test/clib/test_ctonedim.cpp index 4518a4b2b..124d2b245 100644 --- a/test/clib/test_ctonedim.cpp +++ b/test/clib/test_ctonedim.cpp @@ -102,7 +102,7 @@ TEST(ctonedim, reacting_surface_from_parts) int gas = thermo_newFromFile("ptcombust.yaml", "gas"); int surf = thermo_newFromFile("ptcombust.yaml", "Pt_surf"); - int kin = kin_newFromFile("ptcombust.yaml", "Pt_surf", surf, gas, -1, -1, -1); + int kin = kin_newFromFile("ptcombust.yaml", "", surf, gas, -1, -1, -1); ASSERT_GE(kin, 0); int ret = reactingsurf_setkineticsmgr(index, kin); @@ -117,7 +117,7 @@ TEST(ctonedim, catcomb_stack) int trans = soln_transport(sol); int surf = thermo_newFromFile("ptcombust.yaml", "Pt_surf"); - int surf_kin = kin_newFromFile("ptcombust.yaml", "Pt_surf", surf, gas, -1, -1, -1); + int surf_kin = kin_newFromFile("ptcombust.yaml", "", surf, gas, -1, -1, -1); // inlet int inlet = inlet_new(); diff --git a/test/clib/test_ctreactor.cpp b/test/clib/test_ctreactor.cpp index 3b6688986..3188c39e1 100644 --- a/test/clib/test_ctreactor.cpp +++ b/test/clib/test_ctreactor.cpp @@ -11,7 +11,7 @@ using namespace Cantera; TEST(ctreactor, reactor_objects) { int thermo = thermo_newFromFile("gri30.yaml", "gri30"); - int kin = kin_newFromFile("gri30.yaml", "gri30", thermo, -1, -1, -1, -1); + int kin = kin_newFromFile("gri30.yaml", "", thermo, -1, -1, -1, -1); int reactor = reactor_new("IdealGasReactor"); ASSERT_GT(reactor, 0); @@ -64,7 +64,7 @@ TEST(ctreactor, reactor_from_parts) string X = "CH4:1.0, O2:2.0, N2:7.52"; int thermo = thermo_newFromFile("gri30.yaml", "gri30"); - int kin = kin_newFromFile("gri30.yaml", "gri30", thermo, -1, -1, -1, -1); + int kin = kin_newFromFile("gri30.yaml", "", thermo, -1, -1, -1, -1); thermo_setMoleFractionsByName(thermo, X.c_str()); thermo_setTemperature(thermo, T); thermo_setPressure(thermo, P); diff --git a/test/general/test_serialization.cpp b/test/general/test_serialization.cpp index 32b91d2dd..dd4a16fad 100644 --- a/test/general/test_serialization.cpp +++ b/test/general/test_serialization.cpp @@ -281,7 +281,7 @@ TEST(YamlWriter, Interface) shared_ptr gas1(newThermo("ptcombust.yaml", "gas")); shared_ptr surf1(newThermo("ptcombust.yaml", "Pt_surf")); vector> phases1{surf1, gas1}; - shared_ptr kin1 = newKinetics(phases1, "ptcombust.yaml", "Pt_surf"); + shared_ptr kin1 = newKinetics(phases1, "ptcombust.yaml"); double T = 900; double P = OneAtm; @@ -301,7 +301,7 @@ TEST(YamlWriter, Interface) shared_ptr gas2(newThermo("generated-ptcombust.yaml", "gas")); shared_ptr surf2(newThermo("generated-ptcombust.yaml", "Pt_surf")); vector> phases2{surf2, gas2}; - shared_ptr kin2 = newKinetics(phases2, "generated-ptcombust.yaml", "Pt_surf"); + shared_ptr kin2 = newKinetics(phases2, "generated-ptcombust.yaml"); auto iface1 = std::dynamic_pointer_cast(surf1); auto iface2 = std::dynamic_pointer_cast(surf2); diff --git a/test/kinetics/kineticsFromScratch.cpp b/test/kinetics/kineticsFromScratch.cpp index 54d1eb9c8..940c802dd 100644 --- a/test/kinetics/kineticsFromScratch.cpp +++ b/test/kinetics/kineticsFromScratch.cpp @@ -27,7 +27,7 @@ public: { vector> th; th.push_back(pp_ref); - kin_ref = newKinetics(th, "../data/kineticsfromscratch.yaml", "ohmech"); + kin_ref = newKinetics(th, "../data/kineticsfromscratch.yaml"); kin.addThermo(pp); kin.init(); @@ -519,7 +519,7 @@ public: , surf(newThermo("sofc.yaml", "metal_surface")) , surf_ref(newThermo("sofc.yaml", "metal_surface")) { - kin_ref = newKinetics({surf_ref, gas_ref}, "sofc.yaml", "metal_surface"); + kin_ref = newKinetics({surf_ref, gas_ref}, "sofc.yaml"); kin.addThermo(surf); kin.addThermo(gas); } @@ -602,7 +602,7 @@ public: p = make_shared(); vector> th; th.push_back(pp_ref); - kin_ref = newKinetics(th, "../data/kineticsfromscratch.yaml", "ohmech"); + kin_ref = newKinetics(th, "../data/kineticsfromscratch.yaml"); kin.addThermo(p); std::vector> S = getSpecies( diff --git a/test/kinetics/rates.cpp b/test/kinetics/rates.cpp index f9e37ae6a..d8b767fcc 100644 --- a/test/kinetics/rates.cpp +++ b/test/kinetics/rates.cpp @@ -16,7 +16,7 @@ public: FracCoeffTest() : therm(newThermo("frac.yaml", "gas")) { - kin = newKinetics({therm}, "frac.yaml", "gas"); + kin = newKinetics({therm}, "frac.yaml"); therm->setState_TPX(2000, 4*OneAtm, "H2O:0.5, OH:.05, H:0.1, O2:0.15, H2:0.2"); kH2O = therm->speciesIndex("H2O"); diff --git a/test_problems/diamondSurf/runDiamond.cpp b/test_problems/diamondSurf/runDiamond.cpp index 0cf2a077d..9b1c8fe43 100644 --- a/test_problems/diamondSurf/runDiamond.cpp +++ b/test_problems/diamondSurf/runDiamond.cpp @@ -28,8 +28,7 @@ int main(int argc, char** argv) size_t nsp_d100 = diamond100->nSpecies(); cout << "Number of species in diamond_100 = " << nsp_d100 << endl; - auto kin = newKinetics({gas, diamond, diamond100}, - "diamond.yaml", "diamond_100"); + auto kin = newKinetics({diamond100, gas, diamond}, "diamond.yaml"); InterfaceKinetics& ikin = dynamic_cast(*kin); size_t nr = kin->nReactions(); cout << "Number of reactions = " << nr << endl; diff --git a/test_problems/surfSolverTest/results2_blessed.txt b/test_problems/surfSolverTest/results2_blessed.txt index 0d95a2810..40306c47c 100644 --- a/test_problems/surfSolverTest/results2_blessed.txt +++ b/test_problems/surfSolverTest/results2_blessed.txt @@ -1,6 +1,6 @@ Gas Temperature = 1.4e+03 Gas Pressure = 1.01e+05 -Gas Phase: gas (0) +Gas Phase: gas (2) Name Conc MoleF SrcRate (kmol/m^3) (kmol/m^2/s) 0 H 8.69602e-05 0.00999001 -2.365619e-03 @@ -15,7 +15,7 @@ Gas Phase: gas (0) 9 O2 8.69602e-06 0.000999001 -9.060640e-07 Sum of gas mole fractions= 1 -Bulk Phase: soot (10) +Bulk Phase: soot (11) Bulk Temperature = 1.4e+03 Bulk Pressure = 1.01e+05 Name Conc MoleF SrcRate @@ -28,14 +28,14 @@ Density of bulk phase = 3.52e+03 kg / m^3 = 3.52 gm / cm^3 Sum of bulk mole fractions= 1 -Surface Phase: soot_interface (11) +Surface Phase: soot_interface (0) Surface Temperature = 1.4e+03 Surface Pressure = 1.01e+05 Name Coverage SrcRate 0 Csoot-* 0.0184677 0.000000e+00 1 Csoot-H 0.981532 0.000000e+00 Sum of coverages = 1 -Surface Phase: soot_interface (11) +Surface Phase: soot_interface (0) Surface Temperature = 1.4e+03 Surface Pressure = 1.01e+05 Name Coverage SrcRate diff --git a/test_problems/surfSolverTest/surfaceSolver2.cpp b/test_problems/surfSolverTest/surfaceSolver2.cpp index c7b47424e..7c681324a 100644 --- a/test_problems/surfSolverTest/surfaceSolver2.cpp +++ b/test_problems/surfSolverTest/surfaceSolver2.cpp @@ -26,14 +26,15 @@ void printGas(ostream& oooo, shared_ptr gasTP, double x[MSSIZE]; double C[MSSIZE]; oooo.precision(3); - string gasPhaseName = "gas"; + string gasPhaseName = "gas"; gasTP->getMoleFractions(x); gasTP->getConcentrations(C); double Temp = gasTP->temperature(); double p = gasTP->pressure(); oooo << "Gas Temperature = " << Temp << endl; oooo << "Gas Pressure = " << p << endl; - size_t kstart = iKin_ptr->kineticsSpeciesIndex(0, 0); + size_t nGas = iKin_ptr->phaseIndex(gasPhaseName); + size_t kstart = iKin_ptr->kineticsSpeciesIndex(0, nGas); oooo << "Gas Phase: " << gasPhaseName << " " << "(" << kstart << ")" << endl; oooo << " Name " @@ -43,7 +44,7 @@ void printGas(ostream& oooo, shared_ptr gasTP, double sum = 0.0; size_t nspGas = gasTP->nSpecies(); for (size_t k = 0; k < nspGas; k++) { - kstart = iKin_ptr->kineticsSpeciesIndex(k, 0); + kstart = iKin_ptr->kineticsSpeciesIndex(k, nGas); fmt::print(oooo, "{:4d} {:>24s} {:14g} {:14g} {:14e}\n", k, gasTP->speciesName(k), C[k], x[k], src[kstart]); sum += x[k]; @@ -63,7 +64,8 @@ void printBulk(ostream& oooo, shared_ptr bulkPhaseTP, string bulkParticlePhaseName = bulkPhaseTP->name(); bulkPhaseTP->getMoleFractions(x); bulkPhaseTP->getConcentrations(C); - size_t kstart = iKin_ptr->kineticsSpeciesIndex(0, 1); + size_t nBulk = iKin_ptr->phaseIndex(bulkParticlePhaseName); + size_t kstart = iKin_ptr->kineticsSpeciesIndex(0, nBulk); double dens = bulkPhaseTP->density(); oooo << "Bulk Phase: " << bulkParticlePhaseName << " " << "(" << kstart << ")" << endl; @@ -80,7 +82,7 @@ void printBulk(ostream& oooo, shared_ptr bulkPhaseTP, const vector_fp& molecW = bulkPhaseTP->molecularWeights(); size_t nspBulk = bulkPhaseTP->nSpecies(); for (size_t k = 0; k < nspBulk; k++) { - kstart = iKin_ptr->kineticsSpeciesIndex(k, 1); + kstart = iKin_ptr->kineticsSpeciesIndex(k, nBulk); fmt::print(oooo, "{:4d} {:>24s} {:14g} {:14g} {:14e}\n", k, bulkPhaseTP->speciesName(k), C[k], x[k], src[kstart]); sum += x[k]; @@ -105,7 +107,8 @@ void printSurf(ostream& oooo, shared_ptr surfPhaseTP, oooo.precision(3); string surfParticlePhaseName = surfPhaseTP->name(); surfPhaseTP->getMoleFractions(x); - size_t kstart = iKin_ptr->kineticsSpeciesIndex(0, 2); + size_t nSurf = iKin_ptr->phaseIndex(surfParticlePhaseName); + size_t kstart = iKin_ptr->kineticsSpeciesIndex(0, nSurf); oooo << "Surface Phase: " << surfParticlePhaseName << " (" << kstart << ")" << endl; double Temp = surfPhaseTP->temperature(); @@ -117,7 +120,7 @@ void printSurf(ostream& oooo, shared_ptr surfPhaseTP, double sum = 0.0; size_t nspSurf = surfPhaseTP->nSpecies(); for (size_t k = 0; k < nspSurf; k++) { - kstart = iKin_ptr->kineticsSpeciesIndex(0, 2); + kstart = iKin_ptr->kineticsSpeciesIndex(k, nSurf); double srcK = src[kstart]; if (fabs(srcK) < 1.0E-8) { srcK = 0.0; @@ -152,8 +155,7 @@ int main(int argc, char** argv) cout << "Number of species in surface phase, " << surfParticlePhaseName << " = " << nsp_d100 << endl; - auto kin = newKinetics({gasTP, bulkPhaseTP, surfPhaseTP}, - infile, surfParticlePhaseName); + auto kin = newKinetics({surfPhaseTP, gasTP, bulkPhaseTP}, infile); auto iKin_ptr = dynamic_pointer_cast(kin); size_t nr = iKin_ptr->nReactions(); cout << "Number of reactions = " << nr << endl; @@ -170,8 +172,7 @@ int main(int argc, char** argv) // create the second InterfaceKinetics object based on the // second surface phase. - auto kin2 = newKinetics({gasTP, bulkPhaseTP, surfPhaseTP2}, - infile, surfParticlePhaseName); + auto kin2 = newKinetics({surfPhaseTP2, gasTP, bulkPhaseTP}, infile); auto iKin2_ptr = dynamic_pointer_cast(kin2); nr = iKin_ptr->nReactions(); cout << "Number of reactions = " << nr << endl;