diff --git a/include/cantera/equil/MultiPhase.h b/include/cantera/equil/MultiPhase.h
index 03cb31fbcf8..904a2b8ce39 100644
--- a/include/cantera/equil/MultiPhase.h
+++ b/include/cantera/equil/MultiPhase.h
@@ -670,8 +670,8 @@ inline std::ostream& operator<<(std::ostream& s, MultiPhase& x)
 {
     x.updatePhases();
     for (size_t ip = 0; ip < x.nPhases(); ip++) {
-        if (x.phase(ip).name() != "") {
-            s << "*************** " << x.phase(ip).name() << " *****************" << std::endl;
+        if (x.phase(ip).id() != "") {
+            s << "*************** " << x.phase(ip).id() << " *****************" << std::endl;
         } else {
             s << "*************** Phase " << ip << " *****************" << std::endl;
         }
diff --git a/include/cantera/equil/vcs_VolPhase.h b/include/cantera/equil/vcs_VolPhase.h
index 5a4945d013e..f7fadc5aadf 100644
--- a/include/cantera/equil/vcs_VolPhase.h
+++ b/include/cantera/equil/vcs_VolPhase.h
@@ -92,11 +92,11 @@ class vcs_VolPhase
     /*!
      * @param phaseNum    index of the phase in the vcs problem
      * @param numSpecies  Number of species in the phase
-     * @param phaseName   String name for the phase
+     * @param phaseID     String identifier for the phase
      * @param molesInert  kmoles of inert in the phase (defaults to zero)
      */
     void resize(const size_t phaseNum, const size_t numSpecies,
-                const size_t numElem, const char* const phaseName,
+                const size_t numElem, const char* const phaseID,
                 const double molesInert = 0.0);
 
     void elemResize(const size_t numElemConstraints);
@@ -623,10 +623,13 @@ class vcs_VolPhase
     size_t m_numSpecies;
 
 public:
-    //! String name for the phase
-    std::string PhaseName;
+    //! String identifier for the phase
+    std::string PhaseID();
 
 private:
+    //! String identifier of the phase
+    std::string m_phaseID;
+
     //!  Total moles of inert in the phase
     double m_totalMolesInert;
 
diff --git a/include/cantera/equil/vcs_solve.h b/include/cantera/equil/vcs_solve.h
index 9727e0efb14..b3bc4d00f20 100644
--- a/include/cantera/equil/vcs_solve.h
+++ b/include/cantera/equil/vcs_solve.h
@@ -527,7 +527,7 @@ class VCS_SOLVE
      *       in the private data structure. All references to the species
      *       properties must employ the ind[] index vector.
      *    4. Initialization of arrays to zero.
-     *    5. Check to see if the problem is well posed (If all the element 
+     *    5. Check to see if the problem is well posed (If all the element
      *       abundances are zero, the algorithm will fail)
      *
      * @param printLvl Print level of the routine
@@ -1354,7 +1354,7 @@ class VCS_SOLVE
     vector_int m_speciesStatus;
 
     //! Mapping from the species number to the phase number
-    std::vector<size_t> m_phaseID;
+    std::vector<size_t> m_phaseNum;
 
     //! Boolean indicating whether a species belongs to a single-species phase
     // vector<bool> can't be used here because it doesn't work with std::swap
diff --git a/include/cantera/thermo/Phase.h b/include/cantera/thermo/Phase.h
index a82e69d7842..dff2c64588f 100644
--- a/include/cantera/thermo/Phase.h
+++ b/include/cantera/thermo/Phase.h
@@ -66,15 +66,11 @@ namespace Cantera
  * operate on a state vector, which is in general of length (2 + nSpecies()).
  * The first two entries of the state vector are temperature and density.
  *
- * A species name may be referred to via three methods:
- *
- *    -   "speciesName"
- *    -   "PhaseId:speciesName"
- *    -   "phaseName:speciesName"
- *    .
- *
- * The first two methods of naming may not yield a unique species within
- * complicated assemblies of %Cantera Phases.
+ * A species name is referred to via "speciesName", which is unique within a
+ * given phase. Note that within multiphase mixtures ("MultiPhase"), both a
+ * phase name/index as well as species name are required to access information
+ * about a species in a particular phase. For surfaces, the species names are
+ * unique among the phases.
  *
  * @todo
  *   - Make the concept of saving state vectors more general, so that it can
diff --git a/interfaces/cython/cantera/test/test_thermo.py b/interfaces/cython/cantera/test/test_thermo.py
index 19621a73e26..4a57f2ca1e4 100644
--- a/interfaces/cython/cantera/test/test_thermo.py
+++ b/interfaces/cython/cantera/test/test_thermo.py
@@ -304,7 +304,6 @@ def test_name(self):
 
         self.phase.name = 'something'
         self.assertEqual(self.phase.name, 'something')
-        self.assertIn('something', self.phase.report())
 
     def test_ID(self):
         self.assertEqual(self.phase.ID, 'ohmech')
@@ -312,6 +311,7 @@ def test_ID(self):
         self.phase.ID = 'something'
         self.assertEqual(self.phase.ID, 'something')
         self.assertEqual(self.phase.name, 'ohmech')
+        self.assertIn('something', self.phase.report())
 
     def test_badLength(self):
         X = np.zeros(5)
diff --git a/src/clib/ct.cpp b/src/clib/ct.cpp
index bc3e0631d82..2d0796b407a 100644
--- a/src/clib/ct.cpp
+++ b/src/clib/ct.cpp
@@ -287,7 +287,7 @@ extern "C" {
     int thermo_getName(int n, size_t lennm, char* nm)
     {
         try {
-            return static_cast<int>(copyString(ThermoCabinet::item(n).name(), nm, lennm));
+            return static_cast<int>(copyString(ThermoCabinet::item(n).id(), nm, lennm));
         } catch (...) {
             return handleAllExceptions(-1, ERR);
         }
@@ -296,7 +296,7 @@ extern "C" {
     int thermo_setName(int n, const char* nm)
     {
         try {
-            ThermoCabinet::item(n).setName(nm);
+            ThermoCabinet::item(n).setID(nm);
             return 0;
         } catch (...) {
             return handleAllExceptions(-1, ERR);
diff --git a/src/equil/MultiPhase.cpp b/src/equil/MultiPhase.cpp
index dc5f9c48201..f31c4404d60 100644
--- a/src/equil/MultiPhase.cpp
+++ b/src/equil/MultiPhase.cpp
@@ -55,7 +55,7 @@ void MultiPhase::addPhase(ThermoPhase* p, doublereal moles)
 
     if (!p->compatibleWithMultiPhase()) {
         throw CanteraError("MultiPhase::addPhase", "Phase '{}'' is not "
-            "compatible with MultiPhase equilibrium solver", p->name());
+            "compatible with MultiPhase equilibrium solver", p->id());
     }
 
     // save the pointer to the phase object
diff --git a/src/equil/MultiPhaseEquil.cpp b/src/equil/MultiPhaseEquil.cpp
index 3e8c7458829..490ca698642 100644
--- a/src/equil/MultiPhaseEquil.cpp
+++ b/src/equil/MultiPhaseEquil.cpp
@@ -706,7 +706,7 @@ void MultiPhaseEquil::reportCSV(const std::string& reportFile)
         ThermoPhase& tref = m_mix->phase(iphase);
         ThermoPhase* tp = &tref;
         tp->getMoleFractions(&mf[istart]);
-        string phaseName = tref.name();
+        string phaseID = tref.id();
         double TMolesPhase = phaseMoles(iphase);
         size_t nSpecies = tref.nSpecies();
         activity.resize(nSpecies, 0.0);
@@ -733,7 +733,7 @@ void MultiPhaseEquil::reportCSV(const std::string& reportFile)
             tp->getChemPotentials(mu.data());
 
             if (iphase == 0) {
-                fprintf(FP,"        Name,      Phase,  PhaseMoles,  Mole_Fract, "
+                fprintf(FP,"        Name,    PhaseID,  PhaseMoles,  Mole_Fract, "
                         "Molalities,  ActCoeff,   Activity,"
                         "ChemPot_SS0,   ChemPot,   mole_num,       PMVol, Phase_Volume\n");
 
@@ -746,7 +746,7 @@ void MultiPhaseEquil::reportCSV(const std::string& reportFile)
                 fprintf(FP,"%12s, %11s, %11.3e, %11.3e, %11.3e, %11.3e, %11.3e,"
                         "%11.3e, %11.3e, %11.3e, %11.3e, %11.3e\n",
                         sName.c_str(),
-                        phaseName.c_str(), TMolesPhase,
+                        phaseID.c_str(), TMolesPhase,
                         mf[istart + k], molalities[k], ac[k], activity[k],
                         mu0[k]*1.0E-6, mu[k]*1.0E-6,
                         mf[istart + k] * TMolesPhase,
@@ -754,7 +754,7 @@ void MultiPhaseEquil::reportCSV(const std::string& reportFile)
             }
         } else {
             if (iphase == 0) {
-                fprintf(FP,"        Name,       Phase,  PhaseMoles,  Mole_Fract,  "
+                fprintf(FP,"        Name,     PhaseID,  PhaseMoles,  Mole_Fract,  "
                         "Molalities,   ActCoeff,    Activity,"
                         "  ChemPotSS0,     ChemPot,   mole_num,       PMVol, Phase_Volume\n");
 
@@ -770,7 +770,7 @@ void MultiPhaseEquil::reportCSV(const std::string& reportFile)
                 fprintf(FP,"%12s, %11s, %11.3e, %11.3e, %11.3e, %11.3e, %11.3e, "
                         "%11.3e, %11.3e,% 11.3e, %11.3e, %11.3e\n",
                         sName.c_str(),
-                        phaseName.c_str(), TMolesPhase,
+                        phaseID.c_str(), TMolesPhase,
                         mf[istart + k], molalities[k], ac[k],
                         activity[k], mu0[k]*1.0E-6, mu[k]*1.0E-6,
                         mf[istart + k] * TMolesPhase,
diff --git a/src/equil/vcs_Gibbs.cpp b/src/equil/vcs_Gibbs.cpp
index d027c173e4d..d9ec1d8255b 100644
--- a/src/equil/vcs_Gibbs.cpp
+++ b/src/equil/vcs_Gibbs.cpp
@@ -42,7 +42,7 @@ double VCS_SOLVE::vcs_GibbsPhase(size_t iphase, const double* const w,
     double g = 0.0;
     double phaseMols = 0.0;
     for (size_t kspec = 0; kspec < m_numSpeciesRdc; ++kspec) {
-        if (m_phaseID[kspec] == iphase && m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) {
+        if (m_phaseNum[kspec] == iphase && m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) {
             g += w[kspec] * fe[kspec];
             phaseMols += w[kspec];
         }
diff --git a/src/equil/vcs_MultiPhaseEquil.cpp b/src/equil/vcs_MultiPhaseEquil.cpp
index b70c732e4b9..4786dcfa684 100644
--- a/src/equil/vcs_MultiPhaseEquil.cpp
+++ b/src/equil/vcs_MultiPhaseEquil.cpp
@@ -469,7 +469,7 @@ int vcs_MultiPhaseEquil::equilibrate_TP(int estimateEquil,
             } else {
                 plogf("  %15.3e   %15.3e  ", m_mix->speciesMoles(i), m_mix->moleFraction(i));
                 if (m_mix->speciesMoles(i) <= 0.0) {
-                    size_t iph = m_vsolve.m_phaseID[i];
+                    size_t iph = m_vsolve.m_phaseNum[i];
                     vcs_VolPhase* VPhase = m_vsolve.m_VolPhaseList[iph].get();
                     if (VPhase->nSpecies() > 1) {
                         plogf("     -1.000e+300\n");
@@ -533,7 +533,7 @@ void vcs_MultiPhaseEquil::reportCSV(const std::string& reportFile)
 
     for (size_t iphase = 0; iphase < nphase; iphase++) {
         ThermoPhase& tref = m_mix->phase(iphase);
-        string phaseName = tref.name();
+        string phaseID = tref.id();
         vcs_VolPhase* volP = m_vsolve.m_VolPhaseList[iphase].get();
         double TMolesPhase = volP->totalMoles();
         size_t nSpecies = tref.nSpecies();
@@ -562,7 +562,7 @@ void vcs_MultiPhaseEquil::reportCSV(const std::string& reportFile)
             tref.getChemPotentials(&mu[0]);
 
             if (iphase == 0) {
-                fprintf(FP,"        Name,      Phase,  PhaseMoles,  Mole_Fract, "
+                fprintf(FP,"        Name,    PhaseID,  PhaseMoles,  Mole_Fract, "
                         "Molalities,  ActCoeff,   Activity,"
                         "ChemPot_SS0,   ChemPot,   mole_num,       PMVol, Phase_Volume\n");
 
@@ -575,7 +575,7 @@ void vcs_MultiPhaseEquil::reportCSV(const std::string& reportFile)
                 fprintf(FP,"%12s, %11s, %11.3e, %11.3e, %11.3e, %11.3e, %11.3e,"
                         "%11.3e, %11.3e, %11.3e, %11.3e, %11.3e\n",
                         sName.c_str(),
-                        phaseName.c_str(), TMolesPhase,
+                        phaseID.c_str(), TMolesPhase,
                         tref.moleFraction(k), molalities[k], ac[k], activity[k],
                         mu0[k]*1.0E-6, mu[k]*1.0E-6,
                         tref.moleFraction(k) * TMolesPhase,
@@ -583,7 +583,7 @@ void vcs_MultiPhaseEquil::reportCSV(const std::string& reportFile)
             }
         } else {
             if (iphase == 0) {
-                fprintf(FP,"        Name,       Phase,  PhaseMoles,  Mole_Fract,  "
+                fprintf(FP,"        Name,     PhaseID,  PhaseMoles,  Mole_Fract,  "
                         "Molalities,   ActCoeff,    Activity,"
                         "  ChemPotSS0,     ChemPot,   mole_num,       PMVol, Phase_Volume\n");
 
@@ -599,7 +599,7 @@ void vcs_MultiPhaseEquil::reportCSV(const std::string& reportFile)
                 fprintf(FP,"%12s, %11s, %11.3e, %11.3e, %11.3e, %11.3e, %11.3e, "
                         "%11.3e, %11.3e,% 11.3e, %11.3e, %11.3e\n",
                         sName.c_str(),
-                        phaseName.c_str(), TMolesPhase,
+                        phaseID.c_str(), TMolesPhase,
                         tref.moleFraction(k), molalities[k], ac[k],
                         activity[k], mu0[k]*1.0E-6, mu[k]*1.0E-6,
                         tref.moleFraction(k) * TMolesPhase,
diff --git a/src/equil/vcs_VolPhase.cpp b/src/equil/vcs_VolPhase.cpp
index d3e4b3bc9e4..4b476aef826 100644
--- a/src/equil/vcs_VolPhase.cpp
+++ b/src/equil/vcs_VolPhase.cpp
@@ -28,6 +28,7 @@ vcs_VolPhase::vcs_VolPhase(VCS_SOLVE* owningSolverObject) :
     m_numElemConstraints(0),
     m_elemGlobalIndex(0),
     m_numSpecies(0),
+    m_phaseID("<phase-id>"),
     m_totalMolesInert(0.0),
     m_isIdealSoln(false),
     m_existence(VCS_PHASE_EXIST_NO),
@@ -58,8 +59,13 @@ vcs_VolPhase::~vcs_VolPhase()
     }
 }
 
+std::string vcs_VolPhase::PhaseID()
+{
+    return m_phaseID;
+}
+
 void vcs_VolPhase::resize(const size_t phaseNum, const size_t nspecies,
-                          const size_t numElem, const char* const phaseName,
+                          const size_t numElem, const char* const phaseID,
                           const double molesInert)
 {
     AssertThrowMsg(nspecies > 0, "vcs_VolPhase::resize", "nspecies Error");
@@ -68,17 +74,17 @@ void vcs_VolPhase::resize(const size_t phaseNum, const size_t nspecies,
     m_phiVarIndex = npos;
 
     if (phaseNum == VP_ID_) {
-        if (strcmp(PhaseName.c_str(), phaseName)) {
+        if (strcmp(m_phaseID.c_str(), phaseID)) {
             throw CanteraError("vcs_VolPhase::resize",
-                               "Strings are different: " + PhaseName + " " +
-                               phaseName + " :unknown situation");
+                               "Strings are different: " + m_phaseID + " " +
+                               phaseID + " :unknown situation");
         }
     } else {
         VP_ID_ = phaseNum;
-        if (!phaseName) {
-            PhaseName = fmt::format("Phase_{}", VP_ID_);
+        if (!phaseID) {
+            m_phaseID = fmt::format("Phase_{}", VP_ID_);
         } else {
-            PhaseName = phaseName;
+            m_phaseID = phaseID;
         }
     }
     if (nspecies > 1) {
@@ -589,7 +595,7 @@ void vcs_VolPhase::setPtrThermoPhase(ThermoPhase* tp_ptr)
         if (m_numSpecies != 0) {
             plogf("Warning Nsp != NVolSpeces: %d %d \n", nsp, m_numSpecies);
         }
-        resize(VP_ID_, nsp, nelem, PhaseName.c_str());
+        resize(VP_ID_, nsp, nelem, m_phaseID.c_str());
     }
     TP_ptr->getMoleFractions(&Xmol_[0]);
     creationMoleNumbers_ = Xmol_;
diff --git a/src/equil/vcs_elem.cpp b/src/equil/vcs_elem.cpp
index dd6278be0ba..cf3cb7bec16 100644
--- a/src/equil/vcs_elem.cpp
+++ b/src/equil/vcs_elem.cpp
@@ -87,7 +87,7 @@ void VCS_SOLVE::vcs_elabPhase(size_t iphase, double* const elemAbundPhase)
     for (size_t j = 0; j < m_nelem; ++j) {
         elemAbundPhase[j] = 0.0;
         for (size_t i = 0; i < m_nsp; ++i) {
-            if (m_speciesUnknownType[i] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE && m_phaseID[i] == iphase) {
+            if (m_speciesUnknownType[i] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE && m_phaseNum[i] == iphase) {
                 elemAbundPhase[j] += m_formulaMatrix(i,j) * m_molNumSpecies_old[i];
             }
         }
diff --git a/src/equil/vcs_inest.cpp b/src/equil/vcs_inest.cpp
index f8843ab178c..d6a76abaa15 100644
--- a/src/equil/vcs_inest.cpp
+++ b/src/equil/vcs_inest.cpp
@@ -83,7 +83,7 @@ void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm,
     }
     for (size_t kspec = 0; kspec < m_numComponents; ++kspec) {
         if (m_speciesUnknownType[kspec] == VCS_SPECIES_TYPE_MOLNUM) {
-            m_tPhaseMoles_new[m_phaseID[kspec]] += m_molNumSpecies_old[kspec];
+            m_tPhaseMoles_new[m_phaseNum[kspec]] += m_molNumSpecies_old[kspec];
         }
     }
     double TMolesMultiphase = 0.0;
@@ -103,7 +103,7 @@ void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm,
     for (size_t kspec = 0; kspec < m_numComponents; ++kspec) {
         if (m_speciesUnknownType[kspec] == VCS_SPECIES_TYPE_MOLNUM) {
             if (! m_SSPhase[kspec]) {
-                size_t iph = m_phaseID[kspec];
+                size_t iph = m_phaseNum[kspec];
                 m_feSpecies_new[kspec] += log(m_molNumSpecies_new[kspec] / m_tPhaseMoles_old[iph]);
             }
         } else {
@@ -140,7 +140,7 @@ void VCS_SOLVE::vcs_inest(double* const aw, double* const sa, double* const sm,
         // the phase exists, it stays. If it doesn't exist in the estimate, it
         // doesn't come into existence here.
         if (! m_SSPhase[kspec]) {
-            size_t iph = m_phaseID[kspec];
+            size_t iph = m_phaseNum[kspec];
             if (m_deltaGRxn_new[irxn] > xtphMax[iph]) {
                 m_deltaGRxn_new[irxn] = 0.8 * xtphMax[iph];
             }
diff --git a/src/equil/vcs_phaseStability.cpp b/src/equil/vcs_phaseStability.cpp
index f064f01dd47..ee033f22972 100644
--- a/src/equil/vcs_phaseStability.cpp
+++ b/src/equil/vcs_phaseStability.cpp
@@ -125,7 +125,7 @@ size_t VCS_SOLVE::vcs_popPhaseID(std::vector<size_t> & phasePopPhaseIDs)
         if (existence > 0) {
             if (m_debug_print_lvl >= 2) {
                 plogf("  ---    %18s %5d           NA       %11.3e\n",
-                      Vphase->PhaseName, existence, m_tPhaseMoles_old[iph]);
+                      Vphase->PhaseID(), existence, m_tPhaseMoles_old[iph]);
             }
         } else {
             if (Vphase->m_singleSpecies) {
@@ -154,7 +154,7 @@ size_t VCS_SOLVE::vcs_popPhaseID(std::vector<size_t> & phasePopPhaseIDs)
 
                 if (m_debug_print_lvl >= 2) {
                     plogf("  ---    %18s %5d %10.3g %10.3g %s\n",
-                          Vphase->PhaseName, existence, Fephase,
+                          Vphase->PhaseID(), existence, Fephase,
                           m_tPhaseMoles_old[iph], anote);
                 }
             } else {
@@ -171,13 +171,13 @@ size_t VCS_SOLVE::vcs_popPhaseID(std::vector<size_t> & phasePopPhaseIDs)
                     }
                     if (m_debug_print_lvl >= 2) {
                         plogf("  ---    %18s %5d  %11.3g %11.3g\n",
-                              Vphase->PhaseName, existence, Fephase,
+                              Vphase->PhaseID(), existence, Fephase,
                               m_tPhaseMoles_old[iph]);
                     }
                 } else {
                     if (m_debug_print_lvl >= 2) {
                         plogf("  ---    %18s %5d   blocked  %11.3g\n",
-                              Vphase->PhaseName,
+                              Vphase->PhaseID(),
                               existence, m_tPhaseMoles_old[iph]);
                     }
                 }
@@ -215,7 +215,7 @@ int VCS_SOLVE::vcs_popPhaseRxnStepSizes(const size_t iphasePop)
                    "called for a phase that exists!");
     if (m_debug_print_lvl >= 2) {
         plogf("  ---  vcs_popPhaseRxnStepSizes() called to pop phase %s %d into existence\n",
-              Vphase->PhaseName, iphasePop);
+              Vphase->PhaseID(), iphasePop);
     }
     // Section for a single-species phase
     if (Vphase->m_singleSpecies) {
@@ -307,7 +307,7 @@ int VCS_SOLVE::vcs_popPhaseRxnStepSizes(const size_t iphasePop)
                 }
             } else {
                 if (m_elType[j] == VCS_ELEM_TYPE_ABSPOS) {
-                    size_t jph = m_phaseID[j];
+                    size_t jph = m_phaseNum[j];
                     if ((jph != iphasePop) && (!m_SSPhase[j])) {
                         double fdeltaJ = fabs(deltaJ);
                         if (m_molNumSpecies_old[j] > 0.0) {
diff --git a/src/equil/vcs_prep.cpp b/src/equil/vcs_prep.cpp
index 3da6c7ab37f..82ea1694bf7 100644
--- a/src/equil/vcs_prep.cpp
+++ b/src/equil/vcs_prep.cpp
@@ -16,7 +16,7 @@ void VCS_SOLVE::vcs_SSPhase()
 {
     vector_int numPhSpecies(m_numPhases, 0);
     for (size_t kspec = 0; kspec < m_nsp; ++kspec) {
-        numPhSpecies[m_phaseID[kspec]]++;
+        numPhSpecies[m_phaseNum[kspec]]++;
     }
 
     // Handle the special case of a single species in a phase that has been
@@ -37,7 +37,7 @@ void VCS_SOLVE::vcs_SSPhase()
     // information concerning the phase ID of species. SSPhase = Boolean
     // indicating whether a species is in a single species phase or not.
     for (size_t kspec = 0; kspec < m_nsp; kspec++) {
-        size_t iph = m_phaseID[kspec];
+        size_t iph = m_phaseNum[kspec];
         vcs_VolPhase* Vphase = m_VolPhaseList[iph].get();
         if (Vphase->m_singleSpecies) {
             m_SSPhase[kspec] = true;
@@ -70,7 +70,7 @@ int VCS_SOLVE::vcs_prep(int printLvl)
     }
 
     for (size_t kspec = 0; kspec < m_nsp; ++kspec) {
-        size_t pID = m_phaseID[kspec];
+        size_t pID = m_phaseNum[kspec];
         size_t spPhIndex = m_speciesLocalPhaseIndex[kspec];
         vcs_VolPhase* vPhase = m_VolPhaseList[pID].get();
         vcs_SpeciesProperties* spProp = vPhase->speciesProperty(spPhIndex);
diff --git a/src/equil/vcs_prob.cpp b/src/equil/vcs_prob.cpp
index 0745390655e..efc0c668a80 100644
--- a/src/equil/vcs_prob.cpp
+++ b/src/equil/vcs_prob.cpp
@@ -41,12 +41,12 @@ void VCS_SOLVE::prob_report(int print_lvl)
         plogf("\t\tPres = %g atm\n", pres_atm);
         plogf("\n");
         plogf("             Phase IDs of species\n");
-        plogf("            species     phaseID        phaseName   ");
+        plogf("            species     phaseNum       phaseID     ");
         plogf(" Initial_Estimated_Moles   Species_Type\n");
         for (size_t i = 0; i < m_nsp; i++) {
-            vcs_VolPhase* Vphase = m_VolPhaseList[m_phaseID[i]].get();
-            plogf("%16s      %5d   %16s", m_mix->speciesName(i), m_phaseID[i],
-                  Vphase->PhaseName);
+            vcs_VolPhase* Vphase = m_VolPhaseList[m_phaseNum[i]].get();
+            plogf("%16s      %5d   %16s", m_mix->speciesName(i), m_phaseNum[i],
+                  Vphase->PhaseID());
             if (m_doEstimateEquil >= 0) {
                 plogf("             %-10.5g", m_molNumSpecies_old[i]);
             } else {
@@ -65,13 +65,13 @@ void VCS_SOLVE::prob_report(int print_lvl)
         // Printout of the Phase structure information
         writeline('-', 80, true, true);
         plogf("             Information about phases\n");
-        plogf("  PhaseName    PhaseNum SingSpec  GasPhase   "
+        plogf("  PhaseID      PhaseNum SingSpec  GasPhase   "
               " EqnState    NumSpec");
         plogf("  TMolesInert      TKmoles\n");
 
         for (size_t iphase = 0; iphase < m_numPhases; iphase++) {
             vcs_VolPhase* Vphase = m_VolPhaseList[iphase].get();
-            plogf("%16s %5d %5d %8d ", Vphase->PhaseName,
+            plogf("%16s %5d %5d %8d ", Vphase->PhaseID(),
                   Vphase->VP_ID_, Vphase->m_singleSpecies, Vphase->m_gasPhase);
             plogf("%16s %8d %16e ", Vphase->eos_name(),
                   Vphase->nSpecies(), Vphase->totalMolesInert());
@@ -101,7 +101,7 @@ void VCS_SOLVE::prob_report(int print_lvl)
                 size_t kglob = Vphase->spGlobalIndexVCS(kindex);
                 plogf("%16s ", m_mix->speciesName(kglob));
                 if (kindex == 0) {
-                    plogf("%16s", Vphase->PhaseName);
+                    plogf("%16s", Vphase->PhaseID());
                 } else {
                     plogf("                ");
                 }
diff --git a/src/equil/vcs_report.cpp b/src/equil/vcs_report.cpp
index 348ee7def3c..db6394eba0b 100644
--- a/src/equil/vcs_report.cpp
+++ b/src/equil/vcs_report.cpp
@@ -84,7 +84,7 @@ int VCS_SOLVE::vcs_report(int iconv)
                 plogf(" Inert Gas Species        ");
             } else {
                 plogf(" Inert Species in phase %16s ",
-                      m_VolPhaseList[i]->PhaseName);
+                      m_VolPhaseList[i]->PhaseID());
             }
             plogf("%14.7E     %14.7E    %12.4E\n", TPhInertMoles[i],
                   TPhInertMoles[i] / m_tPhaseMoles_old[i], 0.0);
@@ -162,7 +162,7 @@ int VCS_SOLVE::vcs_report(int iconv)
         plogf(" %10.10s", m_elementName[j]);
     }
     plogf(" |                     |\n");
-    plogf("    PhaseName     |KMolTarget |");
+    plogf("    PhaseID       |KMolTarget |");
     for (size_t j = 0; j < m_nelem; j++) {
         plogf(" %10.3g", m_elemAbundancesGoal[j]);
     }
@@ -171,7 +171,7 @@ int VCS_SOLVE::vcs_report(int iconv)
     for (size_t iphase = 0; iphase < m_numPhases; iphase++) {
         plogf(" %3d ", iphase);
         vcs_VolPhase* VPhase = m_VolPhaseList[iphase].get();
-        plogf("%-12.12s |",VPhase->PhaseName);
+        plogf("%-12.12s |", VPhase->PhaseID());
         plogf("%10.3e |", m_tPhaseMoles_old[iphase]);
         totalMoles += m_tPhaseMoles_old[iphase];
         if (m_tPhaseMoles_old[iphase] != VPhase->totalMoles() &&
@@ -230,7 +230,7 @@ int VCS_SOLVE::vcs_report(int iconv)
     writeline('-', 147, true, true);
     for (size_t i = 0; i < m_nsp; ++i) {
         size_t j = x_order[i].second;
-        size_t pid = m_phaseID[j];
+        size_t pid = m_phaseNum[j];
         plogf(" %-12.12s", m_speciesName[j]);
         plogf(" %14.7E ", m_molNumSpecies_old[j]);
         plogf("%14.7E  ", m_SSfeSpecies[j]);
diff --git a/src/equil/vcs_rxnadj.cpp b/src/equil/vcs_rxnadj.cpp
index 66328cb1ea7..ca441b36a01 100644
--- a/src/equil/vcs_rxnadj.cpp
+++ b/src/equil/vcs_rxnadj.cpp
@@ -60,7 +60,7 @@ size_t VCS_SOLVE::vcs_RxnStepSizes(int& forceComponentCalc, size_t& kSpecial)
                 if (m_deltaGRxn_new[irxn] < -1.0e-4) {
                     // First decide if this species is part of a multiphase that
                     // is nontrivial in size.
-                    size_t iph = m_phaseID[kspec];
+                    size_t iph = m_phaseNum[kspec];
                     double tphmoles = m_tPhaseMoles_old[iph];
                     double trphmoles = tphmoles / m_totalMolNum;
                     vcs_VolPhase* Vphase = m_VolPhaseList[iph].get();
@@ -250,7 +250,7 @@ size_t VCS_SOLVE::vcs_RxnStepSizes(int& forceComponentCalc, size_t& kSpecial)
                             m_deltaMolNumSpecies[j] = dss * m_stoichCoeffRxnMatrix(j,irxn);
                         }
 
-                        iphDel = m_phaseID[k];
+                        iphDel = m_phaseNum[k];
                         kSpecial = k;
 
                         if (k != kspec) {
@@ -315,7 +315,7 @@ double VCS_SOLVE::vcs_Hessian_diag_adj(size_t irxn, double hessianDiag_Ideal)
 double VCS_SOLVE::vcs_Hessian_actCoeff_diag(size_t irxn)
 {
     size_t kspec = m_indexRxnToSpecies[irxn];
-    size_t kph = m_phaseID[kspec];
+    size_t kph = m_phaseNum[kspec];
     double np_kspec = std::max(m_tPhaseMoles_old[kph], 1e-13);
     double* sc_irxn = m_stoichCoeffRxnMatrix.ptrColumn(irxn);
 
@@ -327,14 +327,14 @@ double VCS_SOLVE::vcs_Hessian_actCoeff_diag(size_t irxn)
     for (size_t j = 0; j < m_numComponents; j++) {
         if (!m_SSPhase[j]) {
             for (size_t k = 0; k < m_numComponents; ++k) {
-                if (m_phaseID[k] == m_phaseID[j]) {
-                    double np = m_tPhaseMoles_old[m_phaseID[k]];
+                if (m_phaseNum[k] == m_phaseNum[j]) {
+                    double np = m_tPhaseMoles_old[m_phaseNum[k]];
                     if (np > 0.0) {
                         s += sc_irxn[k] * sc_irxn[j] * m_np_dLnActCoeffdMolNum(j,k) / np;
                     }
                 }
             }
-            if (kph == m_phaseID[j]) {
+            if (kph == m_phaseNum[j]) {
                 s += sc_irxn[j] * (m_np_dLnActCoeffdMolNum(j,kspec) + m_np_dLnActCoeffdMolNum(kspec,j)) / np_kspec;
             }
         }
diff --git a/src/equil/vcs_solve.cpp b/src/equil/vcs_solve.cpp
index 6a05c14a596..2f289a97e1a 100644
--- a/src/equil/vcs_solve.cpp
+++ b/src/equil/vcs_solve.cpp
@@ -104,7 +104,7 @@ VCS_SOLVE::VCS_SOLVE(MultiPhase* mphase, int printLvl) :
     m_speciesStatus.resize(m_nsp, VCS_SPECIES_MAJOR);
 
     m_SSPhase.resize(2*m_nsp, 0);
-    m_phaseID.resize(m_nsp, 0);
+    m_phaseNum.resize(m_nsp, 0);
     m_speciesName.resize(m_nsp);
 
     // space for activity coefficients for all species. Set it equal to one.
@@ -154,7 +154,7 @@ VCS_SOLVE::VCS_SOLVE(MultiPhase* mphase, int printLvl) :
         // Find out the number of species in the phase
         size_t nSpPhase = tPhase->nSpecies();
         // Find out the name of the phase
-        string phaseName = tPhase->name();
+        string phaseID = tPhase->id();
 
         // Call the basic vcs_VolPhase creation routine.
         // Properties set here:
@@ -162,9 +162,9 @@ VCS_SOLVE::VCS_SOLVE(MultiPhase* mphase, int printLvl) :
         //    ->GasPhase = Boolean indicating whether it is a gas phase
         //    ->NumSpecies = number of species in the phase
         //    ->TMolesInert = Inerts in the phase = 0.0 for cantera
-        //    ->PhaseName  = Name of the phase
+        //    ->PhaseID  = String identifier of the phase
         vcs_VolPhase* VolPhase = m_VolPhaseList[iphase].get();
-        VolPhase->resize(iphase, nSpPhase, nelem, phaseName.c_str(), 0.0);
+        VolPhase->resize(iphase, nSpPhase, nelem, phaseID.c_str(), 0.0);
         VolPhase->m_gasPhase = gasPhase;
 
         // Tell the vcs_VolPhase pointer about cantera
@@ -221,7 +221,7 @@ VCS_SOLVE::VCS_SOLVE(MultiPhase* mphase, int printLvl) :
             m_chargeSpecies[kT] = tPhase->charge(k);
 
             // Set the phaseid of the species
-            m_phaseID[kT] = iphase;
+            m_phaseNum[kT] = iphase;
 
             // Transfer the type of unknown
             m_speciesUnknownType[kT] = VolPhase->speciesUnknownType(k);
@@ -342,14 +342,14 @@ VCS_SOLVE::VCS_SOLVE(MultiPhase* mphase, int printLvl) :
         writeline('=', 20);
         writeline('=', 80);
         plogf("             Phase IDs of species\n");
-        plogf("            species     phaseID        phaseName   ");
+        plogf("            species     phaseNum       phaseID     ");
         plogf(" Initial_Estimated_kMols\n");
         for (size_t i = 0; i < m_nsp; i++) {
-            size_t iphase = m_phaseID[i];
+            size_t iphase = m_phaseNum[i];
 
             vcs_VolPhase* VolPhase = m_VolPhaseList[iphase].get();
             plogf("%16s      %5d   %16s", mphase->speciesName(i).c_str(), iphase,
-                  VolPhase->PhaseName.c_str());
+                  VolPhase->PhaseID().c_str());
             if (m_speciesUnknownType[i] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) {
                 plogf("     Volts = %-10.5g\n", m_molNumSpecies_old[i]);
             } else {
@@ -360,12 +360,12 @@ VCS_SOLVE::VCS_SOLVE(MultiPhase* mphase, int printLvl) :
         // Printout of the Phase structure information
         writeline('-', 80, true, true);
         plogf("             Information about phases\n");
-        plogf("  PhaseName    PhaseNum SingSpec GasPhase EqnState NumSpec");
+        plogf("  PhaseID      PhaseNum SingSpec GasPhase EqnState NumSpec");
         plogf("  TMolesInert       Tmoles(kmol)\n");
 
         for (size_t iphase = 0; iphase < m_numPhases; iphase++) {
             vcs_VolPhase* VolPhase = m_VolPhaseList[iphase].get();
-            plogf("%16s %5d %5d %8d %16s %8d %16e ", VolPhase->PhaseName.c_str(),
+            plogf("%16s %5d %5d %8d %16s %8d %16e ", VolPhase->PhaseID().c_str(),
                   VolPhase->VP_ID_, VolPhase->m_singleSpecies,
                   VolPhase->m_gasPhase, VolPhase->eos_name(),
                   VolPhase->nSpecies(), VolPhase->totalMolesInert());
@@ -401,14 +401,14 @@ VCS_SOLVE::VCS_SOLVE(MultiPhase* mphase, int printLvl) :
     // time.
     std::vector<size_t> numPhSp(m_numPhases, 0);
     for (size_t kspec = 0; kspec < m_nsp; kspec++) {
-        size_t iph = m_phaseID[kspec];
+        size_t iph = m_phaseNum[kspec];
         if (iph >= m_numPhases) {
             throw CanteraError("VCS_SOLVE::VCS_SOLVE",
                 "Species to Phase Mapping, PhaseID, has a bad value\n"
-                "\tm_phaseID[{}] = {}\n"
+                "\tm_phaseNum[{}] = {}\n"
                 "Allowed values: 0 to {}", kspec, iph, m_numPhases - 1);
         }
-        m_phaseID[kspec] = m_phaseID[kspec];
+        m_phaseNum[kspec] = m_phaseNum[kspec];
         m_speciesLocalPhaseIndex[kspec] = numPhSp[iph];
         numPhSp[iph]++;
     }
@@ -417,7 +417,7 @@ VCS_SOLVE::VCS_SOLVE(MultiPhase* mphase, int printLvl) :
         if (numPhSp[iph] != Vphase->nSpecies()) {
             throw CanteraError("VCS_SOLVE::VCS_SOLVE",
                 "Number of species in phase {}, {}, doesn't match ({} != {}) [vphase = {}]",
-                ser, iph, Vphase->PhaseName, numPhSp[iph], Vphase->nSpecies(), (size_t) Vphase);
+                ser, iph, Vphase->PhaseID(), numPhSp[iph], Vphase->nSpecies(), (size_t) Vphase);
         }
     }
 
@@ -576,13 +576,13 @@ void VCS_SOLVE::vcs_prob_specifyFully()
         writeline('=', 80);
         plogf("\n");
         plogf("             Phase IDs of species\n");
-        plogf("            species     phaseID        phaseName   ");
+        plogf("            species     phaseNum       phaseID     ");
         plogf(" Initial_Estimated_kMols\n");
         for (size_t i = 0; i < m_nsp; i++) {
-            size_t iphase = m_phaseID[i];
+            size_t iphase = m_phaseNum[i];
             vcs_VolPhase* VolPhase = m_VolPhaseList[iphase].get();
             plogf("%16s      %5d   %16s", m_speciesName[i].c_str(), iphase,
-                  VolPhase->PhaseName.c_str());
+                  VolPhase->PhaseID().c_str());
             if (m_speciesUnknownType[i] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) {
                 plogf("     Volts = %-10.5g\n", m_molNumSpecies_old[i]);
             } else {
@@ -593,12 +593,12 @@ void VCS_SOLVE::vcs_prob_specifyFully()
         // Printout of the Phase structure information
         writeline('-', 80, true, true);
         plogf("             Information about phases\n");
-        plogf("  PhaseName    PhaseNum SingSpec GasPhase EqnState NumSpec");
+        plogf("  PhaseID      PhaseNum SingSpec GasPhase EqnState NumSpec");
         plogf("  TMolesInert       Tmoles(kmol)\n");
 
         for (size_t iphase = 0; iphase < m_numPhases; iphase++) {
             vcs_VolPhase* VolPhase = m_VolPhaseList[iphase].get();
-            plogf("%16s %5d %5d %8d %16s %8d %16e ", VolPhase->PhaseName.c_str(),
+            plogf("%16s %5d %5d %8d %16s %8d %16e ", VolPhase->PhaseID().c_str(),
                   VolPhase->VP_ID_, VolPhase->m_singleSpecies,
                   VolPhase->m_gasPhase, VolPhase->eos_name(),
                   VolPhase->nSpecies(), VolPhase->totalMolesInert());
diff --git a/src/equil/vcs_solve_TP.cpp b/src/equil/vcs_solve_TP.cpp
index f5333ab3bcd..5d1bbd14e60 100644
--- a/src/equil/vcs_solve_TP.cpp
+++ b/src/equil/vcs_solve_TP.cpp
@@ -30,7 +30,7 @@ void VCS_SOLVE::checkDelta1(double* const dsLocal,
     vector_fp dchange(m_numPhases, 0.0);
     for (size_t k = 0; k < kspec; k++) {
         if (m_speciesUnknownType[k] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) {
-            size_t iph = m_phaseID[k];
+            size_t iph = m_phaseNum[k];
             dchange[iph] += dsLocal[k];
         }
     }
@@ -124,7 +124,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit)
             for (size_t j = 0; j < m_nelem; ++j) {
                 plogf("% -7.3g ", m_formulaMatrix(i,j));
             }
-            plogf("  %3d  ", m_phaseID[i]);
+            plogf("  %3d  ", m_phaseNum[i]);
             writeline(' ', std::max(55-int(m_nelem)*8, 0), false);
             plogf("%12.5E  %12.5E", RT * m_SSfeSpecies[i], m_molNumSpecies_old[i]);
             if (m_speciesUnknownType[i] == VCS_SPECIES_TYPE_MOLNUM) {
@@ -141,7 +141,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit)
         if (m_molNumSpecies_old[i] < 0.0) {
             plogf("On Input species %-12s has a negative MF, setting it small\n",
                   m_speciesName[i]);
-            size_t iph = m_phaseID[i];
+            size_t iph = m_phaseNum[i];
             double tmp = m_tPhaseMoles_old[iph] * VCS_RELDELETE_SPECIES_CUTOFF * 10;
             tmp = std::max(tmp, VCS_DELETE_MINORSPECIES_CUTOFF*10.);
             m_molNumSpecies_old[i] = tmp;
@@ -254,7 +254,7 @@ int VCS_SOLVE::vcs_solve_TP(int print_lvl, int printDetails, int maxit)
         if (m_SSPhase[kspec]) {
             m_molNumSpecies_new[kspec] = 1.0;
         } else {
-            size_t iph = m_phaseID[kspec];
+            size_t iph = m_phaseNum[kspec];
             if (m_tPhaseMoles_old[iph] != 0.0) {
                 m_molNumSpecies_new[kspec] = m_molNumSpecies_old[kspec] / m_tPhaseMoles_old[iph];
             } else {
@@ -426,7 +426,7 @@ void VCS_SOLVE::solve_tp_inner(size_t& iti, size_t& it1,
         debuglog("   --- Main Loop Treatment -> Circumvented due to Phase Deletion\n", m_debug_print_lvl >= 2);
         for (size_t k = 0; k < m_nsp; k++) {
             m_molNumSpecies_new[k] = m_molNumSpecies_old[k] + m_deltaMolNumSpecies[k];
-            size_t iph = m_phaseID[k];
+            size_t iph = m_phaseNum[k];
             m_tPhaseMoles_new[iph] += m_deltaMolNumSpecies[k];
         }
         if (kspec >= m_numComponents) {
@@ -466,7 +466,7 @@ void VCS_SOLVE::solve_tp_inner(size_t& iti, size_t& it1,
         for (size_t irxn = 0; irxn < m_numRxnRdc; irxn++) {
             size_t kspec = m_indexRxnToSpecies[irxn];
             double* sc_irxn = m_stoichCoeffRxnMatrix.ptrColumn(irxn);
-            size_t iph = m_phaseID[kspec];
+            size_t iph = m_phaseNum[kspec];
             vcs_VolPhase* Vphase = m_VolPhaseList[iph].get();
             ANOTE[0] = '\0';
             double dx;
@@ -710,7 +710,7 @@ void VCS_SOLVE::solve_tp_inner(size_t& iti, size_t& it1,
                         } else {
                             // We are going to zero the single species phase.
                             // Set the existence flag
-                            iph = m_phaseID[kspec];
+                            iph = m_phaseNum[kspec];
                             Vphase = m_VolPhaseList[iph].get();
                             sprintf(ANOTE, "zeroing out SS phase: ");
 
@@ -820,7 +820,7 @@ void VCS_SOLVE::solve_tp_inner(size_t& iti, size_t& it1,
             } else if (m_deltaMolNumSpecies[k] < 0.0) {
                 // If we are here, we then do a step which violates element
                 // conservation.
-                size_t iph = m_phaseID[k];
+                size_t iph = m_phaseNum[k];
                 m_deltaPhaseMoles[iph] -= m_deltaMolNumSpecies[k];
                 m_deltaMolNumSpecies[k] = 0.0;
             }
@@ -974,7 +974,7 @@ void VCS_SOLVE::solve_tp_inner(size_t& iti, size_t& it1,
         writeline('-', 50);
         for (size_t iph = 0; iph < m_numPhases; iph++) {
             vcs_VolPhase* Vphase = m_VolPhaseList[iph].get();
-            plogf("   ---   %18s = %15.7E\n", Vphase->PhaseName, m_tPhaseMoles_new[iph]);
+            plogf("   ---   %18s = %15.7E\n", Vphase->PhaseID(), m_tPhaseMoles_new[iph]);
         }
         plogf("   ");
         writeline('-', 103);
@@ -1355,7 +1355,7 @@ double VCS_SOLVE::vcs_minor_alt_calc(size_t kspec, size_t irxn, bool* do_delete,
 {
     double w_kspec = m_molNumSpecies_old[kspec];
     double dg_irxn = m_deltaGRxn_old[irxn];
-    size_t iph = m_phaseID[kspec];
+    size_t iph = m_phaseNum[kspec];
 
     *do_delete = false;
     if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) {
@@ -1464,14 +1464,14 @@ int VCS_SOLVE::delta_species(const size_t kspec, double* const delta_ptr)
         // ok, we found a positive dx. implement it.
         *delta_ptr = dx;
         m_molNumSpecies_old[kspec] += dx;
-        size_t iph = m_phaseID[kspec];
+        size_t iph = m_phaseNum[kspec];
         m_tPhaseMoles_old[iph] += dx;
         vcs_setFlagsVolPhase(iph, false, VCS_STATECALC_OLD);
 
         for (size_t j = 0; j < m_numComponents; ++j) {
             double tmp = sc_irxn[j] * dx;
             if (tmp != 0.0) {
-                iph = m_phaseID[j];
+                iph = m_phaseNum[j];
                 m_molNumSpecies_old[j] += tmp;
                 m_tPhaseMoles_old[iph] += tmp;
                 vcs_setFlagsVolPhase(iph, false, VCS_STATECALC_OLD);
@@ -1504,7 +1504,7 @@ int VCS_SOLVE::vcs_zero_species(const size_t kspec)
 int VCS_SOLVE::vcs_delete_species(const size_t kspec)
 {
     const size_t klast = m_numSpeciesRdc - 1;
-    const size_t iph = m_phaseID[kspec];
+    const size_t iph = m_phaseNum[kspec];
     vcs_VolPhase* const Vphase = m_VolPhaseList[iph].get();
     const size_t irxn = kspec - m_numComponents;
 
@@ -1547,7 +1547,7 @@ int VCS_SOLVE::vcs_delete_species(const size_t kspec)
         bool stillExists = false;
         for (size_t k = 0; k < m_numSpeciesRdc; k++) {
             if (m_speciesUnknownType[k] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE &&
-                m_phaseID[k] == iph && m_molNumSpecies_old[k] > 0.0) {
+                m_phaseNum[k] == iph && m_molNumSpecies_old[k] > 0.0) {
                 stillExists = true;
                 break;
             }
@@ -1564,7 +1564,7 @@ int VCS_SOLVE::vcs_delete_species(const size_t kspec)
 
 void VCS_SOLVE::vcs_reinsert_deleted(size_t kspec)
 {
-    size_t iph = m_phaseID[kspec];
+    size_t iph = m_phaseNum[kspec];
     if (m_debug_print_lvl >= 2) {
         plogf("   --- Add back a deleted species: %-12s\n", m_speciesName[kspec]);
     }
@@ -1594,7 +1594,7 @@ void VCS_SOLVE::vcs_reinsert_deleted(size_t kspec)
         if (Vphase->exists() == VCS_PHASE_EXIST_NO) {
             Vphase->setExistence(VCS_PHASE_EXIST_YES);
             for (size_t k = 0; k < m_nsp; k++) {
-                if (m_phaseID[k] == iph && m_speciesStatus[k] != VCS_SPECIES_DELETED) {
+                if (m_phaseNum[k] == iph && m_speciesStatus[k] != VCS_SPECIES_DELETED) {
                     m_speciesStatus[k] = VCS_SPECIES_MINOR;
                 }
             }
@@ -1621,12 +1621,12 @@ bool VCS_SOLVE::vcs_delete_multiphase(const size_t iph)
     // set the phase existence flag to dead
     Vphase->setTotalMoles(0.0);
     if (m_debug_print_lvl >= 2) {
-        plogf("   --- delete_multiphase %d, %s\n", iph, Vphase->PhaseName);
+        plogf("   --- delete_multiphase %d, %s\n", iph, Vphase->PhaseID());
     }
 
     // Loop over all of the species in the phase.
     for (size_t kspec = m_numComponents; kspec < m_numSpeciesRdc; ++kspec) {
-        if (m_phaseID[kspec] == iph) {
+        if (m_phaseNum[kspec] == iph) {
             if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) {
                 // calculate an extent of rxn, dx, that zeroes out the species.
                 double dx = - m_molNumSpecies_old[kspec];
@@ -1637,7 +1637,7 @@ bool VCS_SOLVE::vcs_delete_multiphase(const size_t iph)
                     successful = false;
                     if (m_debug_print_lvl >= 2) {
                         plogf("   --- delete_multiphase %d, %s ERROR problems deleting species %s\n",
-                              iph, Vphase->PhaseName, m_speciesName[kspec]);
+                              iph, Vphase->PhaseID(), m_speciesName[kspec]);
                         plogf("   ---     delta  attempted: %g  achieved: %g   "
                               "  Zeroing it manually\n", dx, dxTent);
                     }
@@ -1662,7 +1662,7 @@ bool VCS_SOLVE::vcs_delete_multiphase(const size_t iph)
 
     double dxPerm = 0.0, dxPerm2 = 0.0;
     for (size_t kcomp = 0; kcomp < m_numComponents; ++kcomp) {
-        if (m_phaseID[kcomp] == iph) {
+        if (m_phaseNum[kcomp] == iph) {
             if (m_debug_print_lvl >= 2) {
                 plogf("   --- delete_multiphase   One of the species is a component %d - %s with mole number %g\n",
                       kcomp, m_speciesName[kcomp], m_molNumSpecies_old[kcomp]);
@@ -1670,7 +1670,7 @@ bool VCS_SOLVE::vcs_delete_multiphase(const size_t iph)
             if (m_molNumSpecies_old[kcomp] != 0.0) {
                 for (size_t kspec = m_numComponents; kspec < m_numSpeciesRdc; ++kspec) {
                     size_t irxn = kspec - m_numComponents;
-                    if (m_phaseID[kspec] != iph) {
+                    if (m_phaseNum[kspec] != iph) {
                         if (m_stoichCoeffRxnMatrix(kcomp,irxn) != 0.0) {
                             double dxWant = -m_molNumSpecies_old[kcomp] / m_stoichCoeffRxnMatrix(kcomp,irxn);
                             if (dxWant + m_molNumSpecies_old[kspec]  < 0.0) {
@@ -1678,7 +1678,7 @@ bool VCS_SOLVE::vcs_delete_multiphase(const size_t iph)
                             }
                             for (size_t jcomp = 0; kcomp < m_numComponents; ++kcomp) {
                                 if (jcomp != kcomp) {
-                                    if (m_phaseID[jcomp] == iph) {
+                                    if (m_phaseNum[jcomp] == iph) {
                                         dxPerm = 0.0;
                                     } else {
                                         double dj = dxWant * m_stoichCoeffRxnMatrix(jcomp,irxn);
@@ -1717,7 +1717,7 @@ bool VCS_SOLVE::vcs_delete_multiphase(const size_t iph)
     // zero mole numbers are still considered active. Whether the phase pops
     // back into existence or not is checked as part of the main iteration loop.
     for (size_t kspec = m_numSpeciesRdc; kspec < m_nsp; ++kspec) {
-        if (m_phaseID[kspec] == iph) {
+        if (m_phaseNum[kspec] == iph) {
             m_molNumSpecies_old[kspec] = 0.0;
             m_molNumSpecies_new[kspec] = 0.0;
             m_deltaMolNumSpecies[kspec] = 0.0;
@@ -1761,7 +1761,7 @@ int VCS_SOLVE::vcs_recheck_deleted()
     // deleted species. Then, calculate Delta G for for formation reactions.
     // Note: fe[] here includes everything except for the ln(x[i]) term
     for (size_t kspec = m_numSpeciesRdc; kspec < m_nsp; ++kspec) {
-        size_t iph = m_phaseID[kspec];
+        size_t iph = m_phaseNum[kspec];
         m_feSpecies_new[kspec] = (m_SSfeSpecies[kspec] + log(m_actCoeffSpecies_old[kspec])
                                   - m_lnMnaughtSpecies[kspec]
                                   + m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iph]);
@@ -1803,7 +1803,7 @@ int VCS_SOLVE::vcs_recheck_deleted()
     int npb = 0;
     for (size_t irxn = m_numRxnRdc; irxn < m_numRxnTot; ++irxn) {
         size_t kspec = m_indexRxnToSpecies[irxn];
-        size_t iph = m_phaseID[kspec];
+        size_t iph = m_phaseNum[kspec];
         if (m_tPhaseMoles_old[iph] == 0.0) {
             if (m_deltaGRxn_new[irxn] < 0.0) {
                 vcs_reinsert_deleted(kspec);
@@ -1835,7 +1835,7 @@ size_t VCS_SOLVE::vcs_add_all_deleted()
     m_molNumSpecies_new = m_molNumSpecies_old;
     for (int cits = 0; cits < 3; cits++) {
         for (size_t kspec = m_numSpeciesRdc; kspec < m_nsp; ++kspec) {
-            size_t iph = m_phaseID[kspec];
+            size_t iph = m_phaseNum[kspec];
             vcs_VolPhase* Vphase = m_VolPhaseList[iph].get();
             if (m_molNumSpecies_new[kspec] == 0.0) {
                 m_molNumSpecies_new[kspec] = VCS_DELETE_MINORSPECIES_CUTOFF * 1.0E-10;
@@ -1852,7 +1852,7 @@ size_t VCS_SOLVE::vcs_add_all_deleted()
         vcs_deltag(0, true, VCS_STATECALC_NEW);
         for (size_t irxn = m_numRxnRdc; irxn < m_numRxnTot; ++irxn) {
             size_t kspec = m_indexRxnToSpecies[irxn];
-            size_t iph = m_phaseID[kspec];
+            size_t iph = m_phaseNum[kspec];
             if (m_tPhaseMoles_old[iph] > 0.0) {
                 double maxDG = std::min(m_deltaGRxn_new[irxn], 690.0);
                 double dx = m_tPhaseMoles_old[iph] * exp(- maxDG);
@@ -1862,7 +1862,7 @@ size_t VCS_SOLVE::vcs_add_all_deleted()
     }
     for (size_t irxn = m_numRxnRdc; irxn < m_numRxnTot; ++irxn) {
         size_t kspec = m_indexRxnToSpecies[irxn];
-        size_t iph = m_phaseID[kspec];
+        size_t iph = m_phaseNum[kspec];
         if (m_tPhaseMoles_old[iph] > 0.0) {
             double dx = m_molNumSpecies_new[kspec];
             size_t retn = delta_species(kspec, &dx);
@@ -1897,7 +1897,7 @@ size_t VCS_SOLVE::vcs_add_all_deleted()
     size_t retn = 0;
     for (size_t irxn = m_numRxnRdc; irxn < m_numRxnTot; ++irxn) {
         size_t kspec = m_indexRxnToSpecies[irxn];
-        size_t iph = m_phaseID[kspec];
+        size_t iph = m_phaseNum[kspec];
         if (m_tPhaseMoles_old[iph] > 0.0 && fabs(m_deltaGRxn_old[irxn]) > m_tolmin) {
             if (((m_molNumSpecies_old[kspec] * exp(-m_deltaGRxn_old[irxn])) >
                     VCS_DELETE_MINORSPECIES_CUTOFF) ||
@@ -2481,11 +2481,11 @@ int VCS_SOLVE::vcs_basopt(const bool doJustComponents, double aw[], double sa[],
     for (size_t irxn = 0; irxn < m_numRxnTot; ++irxn) {
         double* scrxn_ptr = m_stoichCoeffRxnMatrix.ptrColumn(irxn);
         size_t kspec = m_indexRxnToSpecies[irxn];
-        size_t iph = m_phaseID[kspec];
+        size_t iph = m_phaseNum[kspec];
         m_deltaMolNumPhase(iph,irxn) = 1.0;
         m_phaseParticipation(iph,irxn)++;
         for (size_t j = 0; j < ncTrial; ++j) {
-            iph = m_phaseID[j];
+            iph = m_phaseNum[j];
             if (fabs(scrxn_ptr[j]) <= 1.0e-6) {
                 scrxn_ptr[j] = 0.0;
             } else {
@@ -2561,7 +2561,7 @@ int VCS_SOLVE::vcs_species_type(const size_t kspec) const
         return VCS_SPECIES_INTERFACIALVOLTAGE;
     }
 
-    size_t iph = m_phaseID[kspec];
+    size_t iph = m_phaseNum[kspec];
     int irxn = int(kspec) - int(m_numComponents);
     int phaseExist = m_VolPhaseList[iph]->exists();
 
@@ -2617,7 +2617,7 @@ int VCS_SOLVE::vcs_species_type(const size_t kspec) const
                             }
                         }
                     } else if (negChangeComp < 0.0) {
-                        if (m_VolPhaseList[m_phaseID[j]]->exists() <= 0) {
+                        if (m_VolPhaseList[m_phaseNum[j]]->exists() <= 0) {
                             if (m_debug_print_lvl >= 2) {
                                 plogf("   ---   %s is prevented from popping into existence because"
                                       " a needed component %s is in a zeroed-phase that would be "
@@ -2759,7 +2759,7 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc,
     }
     for (size_t kspec = 0; kspec < m_nsp; kspec++) {
         if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) {
-            size_t iph = m_phaseID[kspec];
+            size_t iph = m_phaseNum[kspec];
             tlogMoles[iph] += molNum[kspec];
         }
     }
@@ -2804,7 +2804,7 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc,
     // with loops over either the major noncomponent species or the minor
     // noncomponent species.
     for (size_t kspec = l1; kspec < l2; ++kspec) {
-        size_t iphase = m_phaseID[kspec];
+        size_t iphase = m_phaseNum[kspec];
         if (m_speciesUnknownType[kspec] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) {
             AssertThrowMsg(molNum[kspec] == m_phasePhi[iphase], "VCS_SOLVE::vcs_dfe",
                 "We have an inconsistency!");
@@ -2822,11 +2822,11 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc,
                                    + m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase];
             } else {
                 if (molNum[kspec] <= VCS_DELETE_MINORSPECIES_CUTOFF) {
-                    size_t iph = m_phaseID[kspec];
+                    size_t iph = m_phaseNum[kspec];
                     if (tPhMoles_ptr[iph] > 0.0) {
                         feSpecies[kspec] = m_SSfeSpecies[kspec]
                                            + log(actCoeff_ptr[kspec] * VCS_DELETE_MINORSPECIES_CUTOFF)
-                                           - tlogMoles[m_phaseID[kspec]] - m_lnMnaughtSpecies[kspec]
+                                           - tlogMoles[m_phaseNum[kspec]] - m_lnMnaughtSpecies[kspec]
                                            + m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase];
                     } else {
                         feSpecies[kspec] = m_SSfeSpecies[kspec] - m_lnMnaughtSpecies[kspec]
@@ -2835,7 +2835,7 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc,
                 } else {
                     feSpecies[kspec] = m_SSfeSpecies[kspec]
                                        + log(actCoeff_ptr[kspec] * molNum[kspec])
-                                       - tlogMoles[m_phaseID[kspec]] - m_lnMnaughtSpecies[kspec]
+                                       - tlogMoles[m_phaseNum[kspec]] - m_lnMnaughtSpecies[kspec]
                                        + m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase];
                 }
             }
@@ -2847,7 +2847,7 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc,
         for (size_t irxn = 0; irxn < m_numRxnRdc; ++irxn) {
             size_t kspec = m_indexRxnToSpecies[irxn];
             if (m_speciesStatus[kspec] != VCS_SPECIES_MINOR) {
-                size_t iphase = m_phaseID[kspec];
+                size_t iphase = m_phaseNum[kspec];
                 if (m_speciesUnknownType[kspec] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) {
                     AssertThrowMsg(molNum[kspec] == m_phasePhi[iphase], "VCS_SOLVE::vcs_dfe",
                                    "We have an inconsistency!");
@@ -2865,11 +2865,11 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc,
                                            + m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase];
                     } else {
                         if (molNum[kspec] <= VCS_DELETE_MINORSPECIES_CUTOFF) {
-                            size_t iph = m_phaseID[kspec];
+                            size_t iph = m_phaseNum[kspec];
                             if (tPhMoles_ptr[iph] > 0.0) {
                                 feSpecies[kspec] = m_SSfeSpecies[kspec]
                                                    + log(actCoeff_ptr[kspec] * VCS_DELETE_MINORSPECIES_CUTOFF)
-                                                   - tlogMoles[m_phaseID[kspec]] - m_lnMnaughtSpecies[kspec]
+                                                   - tlogMoles[m_phaseNum[kspec]] - m_lnMnaughtSpecies[kspec]
                                                    + m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase];
                             } else {
                                 feSpecies[kspec] = m_SSfeSpecies[kspec] - m_lnMnaughtSpecies[kspec]
@@ -2878,7 +2878,7 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc,
                         } else {
                             feSpecies[kspec] = m_SSfeSpecies[kspec]
                                                + log(actCoeff_ptr[kspec] * molNum[kspec])
-                                               - tlogMoles[m_phaseID[kspec]] - m_lnMnaughtSpecies[kspec]
+                                               - tlogMoles[m_phaseNum[kspec]] - m_lnMnaughtSpecies[kspec]
                                                + m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase];
                         }
                     }
@@ -2890,7 +2890,7 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc,
         for (size_t irxn = 0; irxn < m_numRxnRdc; ++irxn) {
             size_t kspec = m_indexRxnToSpecies[irxn];
             if (m_speciesStatus[kspec] == VCS_SPECIES_MINOR) {
-                size_t iphase = m_phaseID[kspec];
+                size_t iphase = m_phaseNum[kspec];
                 if (m_speciesUnknownType[kspec] == VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) {
                     AssertThrowMsg(molNum[kspec] == m_phasePhi[iphase], "VCS_SOLVE::vcs_dfe",
                         "We have an inconsistency!");
@@ -2908,11 +2908,11 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc,
                                            + m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase];
                     } else {
                         if (molNum[kspec] <= VCS_DELETE_MINORSPECIES_CUTOFF) {
-                            size_t iph = m_phaseID[kspec];
+                            size_t iph = m_phaseNum[kspec];
                             if (tPhMoles_ptr[iph] > 0.0) {
                                 feSpecies[kspec] = m_SSfeSpecies[kspec]
                                                    + log(actCoeff_ptr[kspec] * VCS_DELETE_MINORSPECIES_CUTOFF)
-                                                   - tlogMoles[m_phaseID[kspec]] - m_lnMnaughtSpecies[kspec]
+                                                   - tlogMoles[m_phaseNum[kspec]] - m_lnMnaughtSpecies[kspec]
                                                    + m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase];
                             } else {
                                 feSpecies[kspec] = m_SSfeSpecies[kspec] - m_lnMnaughtSpecies[kspec]
@@ -2921,7 +2921,7 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc,
                         } else {
                             feSpecies[kspec] = m_SSfeSpecies[kspec]
                                                + log(actCoeff_ptr[kspec] * molNum[kspec])
-                                               - tlogMoles[m_phaseID[kspec]] - m_lnMnaughtSpecies[kspec]
+                                               - tlogMoles[m_phaseNum[kspec]] - m_lnMnaughtSpecies[kspec]
                                                + m_chargeSpecies[kspec] * m_Faraday_dim * m_phasePhi[iphase];
                         }
                     }
@@ -2956,7 +2956,7 @@ double VCS_SOLVE::vcs_tmoles()
     }
     for (size_t i = 0; i < m_nsp; i++) {
         if (m_speciesUnknownType[i] == VCS_SPECIES_TYPE_MOLNUM) {
-            m_tPhaseMoles_old[m_phaseID[i]] += m_molNumSpecies_old[i];
+            m_tPhaseMoles_old[m_phaseNum[i]] += m_molNumSpecies_old[i];
         }
     }
     double sum = 0.0;
@@ -2980,7 +2980,7 @@ void VCS_SOLVE::check_tmoles() const
         double m_tPhaseMoles_old_a = TPhInertMoles[i];
 
         for (size_t k = 0; k < m_nsp; k++) {
-            if (m_speciesUnknownType[k] == VCS_SPECIES_TYPE_MOLNUM && m_phaseID[k] == i) {
+            if (m_speciesUnknownType[k] == VCS_SPECIES_TYPE_MOLNUM && m_phaseNum[k] == i) {
                 m_tPhaseMoles_old_a += m_molNumSpecies_old[k];
             }
         }
@@ -3322,7 +3322,7 @@ void VCS_SOLVE::vcs_printDeltaG(const int stateCalc)
             irxn = kspec - m_numComponents;
         }
         double mfValue = 1.0;
-        size_t iphase = m_phaseID[kspec];
+        size_t iphase = m_phaseNum[kspec];
         const vcs_VolPhase* Vphase = m_VolPhaseList[iphase].get();
         if ((m_speciesStatus[kspec] == VCS_SPECIES_ZEROEDMS) ||
                 (m_speciesStatus[kspec] == VCS_SPECIES_ZEROEDPHASE) ||
@@ -3398,8 +3398,8 @@ void VCS_SOLVE::vcs_switch_pos(const bool ifunc, const size_t k1, const size_t k
     }
 
     // Handle the index pointer in the phase structures first
-    vcs_VolPhase* pv1 = m_VolPhaseList[m_phaseID[k1]].get();
-    vcs_VolPhase* pv2 = m_VolPhaseList[m_phaseID[k2]].get();
+    vcs_VolPhase* pv1 = m_VolPhaseList[m_phaseNum[k1]].get();
+    vcs_VolPhase* pv2 = m_VolPhaseList[m_phaseNum[k2]].get();
     size_t kp1 = m_speciesLocalPhaseIndex[k1];
     size_t kp2 = m_speciesLocalPhaseIndex[k2];
     AssertThrowMsg(pv1->spGlobalIndexVCS(kp1) == k1, "VCS_SOLVE::vcs_switch_pos",
@@ -3418,7 +3418,7 @@ void VCS_SOLVE::vcs_switch_pos(const bool ifunc, const size_t k1, const size_t k
     std::swap(m_feSpecies_old[k1], m_feSpecies_old[k2]);
     std::swap(m_feSpecies_new[k1], m_feSpecies_new[k2]);
     std::swap(m_SSPhase[k1], m_SSPhase[k2]);
-    std::swap(m_phaseID[k1], m_phaseID[k2]);
+    std::swap(m_phaseNum[k1], m_phaseNum[k2]);
     std::swap(m_speciesMapIndex[k1], m_speciesMapIndex[k2]);
     std::swap(m_speciesLocalPhaseIndex[k1], m_speciesLocalPhaseIndex[k2]);
     std::swap(m_actConventionSpecies[k1], m_actConventionSpecies[k2]);
diff --git a/src/thermo/BinarySolutionTabulatedThermo.cpp b/src/thermo/BinarySolutionTabulatedThermo.cpp
index c43a17b67da..1061e9bbc90 100644
--- a/src/thermo/BinarySolutionTabulatedThermo.cpp
+++ b/src/thermo/BinarySolutionTabulatedThermo.cpp
@@ -91,7 +91,7 @@ void BinarySolutionTabulatedThermo::initThermo()
             throw InputFileError("BinarySolutionTabulatedThermo::initThermo",
                 m_input["tabulated-species"],
                 "Species '{}' is not in phase '{}'",
-                m_input["tabulated-species"].asString(), name());
+                m_input["tabulated-species"].asString(), id());
         }
         const AnyMap& table = m_input["tabulated-thermo"].as<AnyMap>();
         vector_fp x = table["mole-fractions"].asVector<double>();
diff --git a/src/thermo/LatticeSolidPhase.cpp b/src/thermo/LatticeSolidPhase.cpp
index 02293b85ba7..484c1e37623 100644
--- a/src/thermo/LatticeSolidPhase.cpp
+++ b/src/thermo/LatticeSolidPhase.cpp
@@ -343,7 +343,7 @@ void LatticeSolidPhase::addLattice(shared_ptr<ThermoPhase> lattice)
 void LatticeSolidPhase::setLatticeStoichiometry(const compositionMap& comp)
 {
     for (size_t i = 0; i < m_lattice.size(); i++) {
-        theta_[i] = getValue(comp, m_lattice[i]->name(), 0.0);
+        theta_[i] = getValue(comp, m_lattice[i]->id(), 0.0);
     }
     // Add in the lattice stoichiometry constraint
     for (size_t i = 1; i < m_lattice.size(); i++) {
diff --git a/src/thermo/MolalityVPSSTP.cpp b/src/thermo/MolalityVPSSTP.cpp
index 1571f1a5bff..d59af294700 100644
--- a/src/thermo/MolalityVPSSTP.cpp
+++ b/src/thermo/MolalityVPSSTP.cpp
@@ -380,8 +380,8 @@ std::string MolalityVPSSTP::report(bool show_thermo, doublereal threshold) const
 {
     fmt::memory_buffer b;
     try {
-        if (name() != "") {
-            format_to(b, "\n  {}:\n", name());
+        if (id() != "") {
+            format_to(b, "\n  {}:\n", id());
         }
         format_to(b, "\n");
         format_to(b, "       temperature    {:12.6g}  K\n", temperature());
diff --git a/src/thermo/Phase.cpp b/src/thermo/Phase.cpp
index 85798d5ec82..20e362864ef 100644
--- a/src/thermo/Phase.cpp
+++ b/src/thermo/Phase.cpp
@@ -206,6 +206,11 @@ size_t Phase::speciesIndex(const std::string& nameStr) const
         std::string pn;
         std::string sn = parseSpeciesName(nameStr, pn);
         if (pn == "" || pn == m_name || pn == m_id) {
+            warn_deprecated("Phase::speciesIndex()",
+                            "Retrieval of species indices via "
+                            "'PhaseId:speciesName' or 'phaseName:"
+                            "speciesName to be removed "
+                            "after Cantera 2.5.");
             it = m_speciesIndices.find(nameStr);
             if (it != m_speciesIndices.end()) {
                 return it->second;
diff --git a/src/thermo/PureFluidPhase.cpp b/src/thermo/PureFluidPhase.cpp
index 2382ac51199..71fafe3a918 100644
--- a/src/thermo/PureFluidPhase.cpp
+++ b/src/thermo/PureFluidPhase.cpp
@@ -385,8 +385,8 @@ void PureFluidPhase::setState_Psat(doublereal p, doublereal x)
 std::string PureFluidPhase::report(bool show_thermo, doublereal threshold) const
 {
     fmt::memory_buffer b;
-    if (name() != "") {
-        format_to(b, "\n  {}:\n", name());
+    if (id() != "") {
+        format_to(b, "\n  {}:\n", id());
     }
     format_to(b, "\n");
     format_to(b, "       temperature    {:12.6g}  K\n", temperature());
diff --git a/src/thermo/ThermoFactory.cpp b/src/thermo/ThermoFactory.cpp
index 387b47ccbe2..239c4af6a3e 100644
--- a/src/thermo/ThermoFactory.cpp
+++ b/src/thermo/ThermoFactory.cpp
@@ -448,6 +448,7 @@ void addSpecies(ThermoPhase& thermo, const AnyValue& names, const AnyValue& spec
 void setupPhase(ThermoPhase& thermo, AnyMap& phaseNode, const AnyMap& rootNode)
 {
     thermo.setName(phaseNode["name"].asString());
+    thermo.setID(phaseNode["name"].asString());
     if (rootNode.hasKey("__file__")) {
         phaseNode["__file__"] = rootNode["__file__"];
     }
diff --git a/src/thermo/ThermoPhase.cpp b/src/thermo/ThermoPhase.cpp
index 8ebbee16eba..044a520493c 100644
--- a/src/thermo/ThermoPhase.cpp
+++ b/src/thermo/ThermoPhase.cpp
@@ -926,8 +926,8 @@ std::string ThermoPhase::report(bool show_thermo, doublereal threshold) const
 {
     fmt::memory_buffer b;
     try {
-        if (name() != "") {
-            format_to(b, "\n  {}:\n", name());
+        if (id() != "") {
+            format_to(b, "\n  {}:\n", id());
         }
         format_to(b, "\n");
         format_to(b, "       temperature    {:12.6g}  K\n", temperature());
diff --git a/test/thermo/phaseConstructors.cpp b/test/thermo/phaseConstructors.cpp
index f9808be6ea4..00f8dc2a34a 100644
--- a/test/thermo/phaseConstructors.cpp
+++ b/test/thermo/phaseConstructors.cpp
@@ -508,14 +508,14 @@ TEST(MargulesVPSSTP, fromScratch)
 TEST(LatticeSolidPhase, fromScratch)
 {
     auto base = make_shared<StoichSubstance>();
-    base->setName("Li7Si3(S)");
+    base->setID("Li7Si3(S)");
     base->setDensity(1390.0);
     auto sLi7Si3 = make_shomate2_species("Li7Si3(S)", "Li:7 Si:3", li7si3_shomate_coeffs);
     base->addSpecies(sLi7Si3);
     base->initThermo();
 
     auto interstital = make_shared<LatticePhase>();
-    interstital->setName("Li7Si3_Interstitial");
+    interstital->setID("Li7Si3_Interstitial");
     auto sLii = make_const_cp_species("Li(i)", "Li:1", 298.15, 0, 2e4, 2e4);
     auto sVac = make_const_cp_species("V(i)", "", 298.15, 8.98e4, 0, 0);
     sLii->extra["molar_volume"] = 0.2;