diff --git a/src/base/ctml.cpp b/src/base/ctml.cpp
index 9ef0197be..db868ffcb 100644
--- a/src/base/ctml.cpp
+++ b/src/base/ctml.cpp
@@ -366,12 +366,12 @@ size_t getFloatArray(const XML_Node& node, vector_fp & v,
         }
         doublereal vv = v.back();
         if (vmin != Undef && vv < vmin - Tiny) {
-            writelog("\nWarning: value "+fp2str(vv)+
-                     " is below lower limit of " +fp2str(vmin)+".\n");
+            writelog("\nWarning: value {} is below lower limit of {}.\n",
+                     vv, vmin);
         }
         if (vmax != Undef && vv > vmax + Tiny) {
-            writelog("\nWarning: value "+fp2str(vv)+
-                     " is above upper limit of " +fp2str(vmin)+".\n");
+            writelog("\nWarning: value {} is above upper limit of {}.\n",
+                     vv, vmax);
         }
     }
     for (size_t n = 0; n < v.size(); n++) {
diff --git a/src/clib/Cabinet.h b/src/clib/Cabinet.h
index c1cd0578c..4a25f5f90 100644
--- a/src/clib/Cabinet.h
+++ b/src/clib/Cabinet.h
@@ -141,7 +141,7 @@ public:
         if (n < data.size()) {
             return *data[n];
         } else {
-            throw Cantera::CanteraError("Cabinet::item","index out of range"+Cantera::int2str(n));
+            throw Cantera::CanteraError("Cabinet::item","index out of range {}", n);
         }
     }
 
diff --git a/src/clib/ctxml.cpp b/src/clib/ctxml.cpp
index 132759c11..b39c1814d 100644
--- a/src/clib/ctxml.cpp
+++ b/src/clib/ctxml.cpp
@@ -287,7 +287,7 @@ extern "C" {
             // array not big enough
             if (n < nv) {
                 throw CanteraError("ctml_getFloatArray",
-                                   "array must be dimensioned at least "+int2str(nv));
+                                   "array must be dimensioned at least {}", nv);
             }
 
             for (size_t i = 0; i < nv; i++) {
diff --git a/src/equil/ChemEquil.cpp b/src/equil/ChemEquil.cpp
index 207610ff5..cdaea65c5 100644
--- a/src/equil/ChemEquil.cpp
+++ b/src/equil/ChemEquil.cpp
@@ -118,12 +118,9 @@ void ChemEquil::initialize(thermo_t& s)
                 // the element should be an electron... if it isn't
                 // print a warning.
                 if (ewt > 1.0e-3) {
-                    writelog(string("WARNING: species "
-                                    +s.speciesName(k)
-                                    +" has "+fp2str(s.nAtoms(k,m))
-                                    +" atoms of element "
-                                    +s.elementName(m)+
-                                    ", but this element is not an electron.\n"));
+                    writelog("WARNING: species {} has {} atoms of element {},"
+                             " but this element is not an electron.\n",
+                             s.speciesName(k), s.nAtoms(k,m), s.elementName(m));
                 }
             }
         }
@@ -174,8 +171,8 @@ void ChemEquil::update(const thermo_t& s)
             m_elementmolefracs[m] += nAtoms(k,m) * m_molefractions[k];
             if (m_molefractions[k] < 0.0) {
                 throw CanteraError("ChemEquil::update",
-                                   "negative mole fraction for "+s.speciesName(k)+
-                                   ": "+fp2str(m_molefractions[k]));
+                                   "negative mole fraction for {}: {}",
+                                   s.speciesName(k), m_molefractions[k]);
             }
         }
         sum += m_elementmolefracs[m];
@@ -403,10 +400,9 @@ int ChemEquil::equilibrate(thermo_t& s, const char* XYstr,
     if (tempFixed) {
         double tfixed = s.temperature();
         if (tfixed > s.maxTemp() + 1.0 || tfixed < s.minTemp() - 1.0) {
-            throw CanteraError("ChemEquil","Specified temperature ("
-                               +fp2str(s.temperature())+" K) outside "
-                               "valid range of "+fp2str(s.minTemp())+" K to "
-                               +fp2str(s.maxTemp())+" K\n");
+            throw CanteraError("ChemEquil::equilibrate", "Specified temperature"
+                               " ({} K) outside valid range of {} K to {} K\n",
+                               s.temperature(), s.minTemp(), s.maxTemp());
         }
     }
 
@@ -670,10 +666,9 @@ int ChemEquil::equilibrate(thermo_t& s, const char* XYstr,
             s.setElementPotentials(m_lambda);
             if (s.temperature() > s.maxTemp() + 1.0 ||
                     s.temperature() < s.minTemp() - 1.0) {
-                writelog("Warning: Temperature ("
-                         +fp2str(s.temperature())+" K) outside "
-                         "valid range of "+fp2str(s.minTemp())+" K to "
-                         +fp2str(s.maxTemp())+" K\n");
+                writelog("Warning: Temperature ({} K) outside valid range of "
+                         "{} K to {} K\n",
+                         s.temperature(), s.minTemp(), s.maxTemp());
             }
             return 0;
         }
@@ -770,9 +765,8 @@ int ChemEquil::equilibrate(thermo_t& s, const char* XYstr,
 
     // no convergence
     s.restoreState(state);
-    throw CanteraError("equilibrate",
-                       "no convergence in "+int2str(options.maxIterations)
-                       +" iterations.");
+    throw CanteraError("ChemEquil::equilibrate",
+                       "no convergence in {} iterations.", options.maxIterations);
 }
 
 
diff --git a/src/equil/MultiPhase.cpp b/src/equil/MultiPhase.cpp
index bf5bb2da7..035f35424 100644
--- a/src/equil/MultiPhase.cpp
+++ b/src/equil/MultiPhase.cpp
@@ -751,7 +751,7 @@ void MultiPhase::equilibrate(const std::string& XY, const std::string& solver,
                                           rtol, max_steps);
             if (ret) {
                 throw CanteraError("MultiPhase::equilibrate",
-                    "VCS solver failed. Return code: " + int2str(ret));
+                    "VCS solver failed. Return code: {}", ret);
             }
             debuglog("VCS solver succeeded\n", log_level);
             return;
diff --git a/src/equil/MultiPhaseEquil.cpp b/src/equil/MultiPhaseEquil.cpp
index 8d8d5d465..f719e54c2 100644
--- a/src/equil/MultiPhaseEquil.cpp
+++ b/src/equil/MultiPhaseEquil.cpp
@@ -179,8 +179,8 @@ doublereal MultiPhaseEquil::equilibrate(int XY, doublereal err,
     }
     if (i >= maxsteps) {
         throw CanteraError("MultiPhaseEquil::equilibrate",
-                           "no convergence in " + int2str(maxsteps) +
-                           " iterations. Error = " + fp2str(error()));
+                           "no convergence in {} iterations. Error = {}",
+                           maxsteps, error());
     }
     finish();
     return error();
diff --git a/src/equil/vcs_MultiPhaseEquil.cpp b/src/equil/vcs_MultiPhaseEquil.cpp
index 107b2562d..484ea8192 100644
--- a/src/equil/vcs_MultiPhaseEquil.cpp
+++ b/src/equil/vcs_MultiPhaseEquil.cpp
@@ -55,7 +55,7 @@ int vcs_MultiPhaseEquil::equilibrate_TV(int XY, doublereal xtarget,
     doublereal Vtarget = m_mix->volume();
     if ((XY != TV) && (XY != HV) && (XY != UV) && (XY != SV)) {
         throw CanteraError("vcs_MultiPhaseEquil::equilibrate_TV",
-                           "Wrong XY flag:" + int2str(XY));
+                           "Wrong XY flag: {}", XY);
     }
     int maxiter = 100;
     int iSuccess = 0;
@@ -159,7 +159,7 @@ int vcs_MultiPhaseEquil::equilibrate_HP(doublereal Htarget,
     int iSuccess;
     if (XY != HP && XY != UP) {
         throw CanteraError("vcs_MultiPhaseEquil::equilibrate_HP",
-                           "Wrong XP" + int2str(XY));
+                           "Wrong XP", XY);
     }
     int strt = estimateEquil;
 
@@ -907,7 +907,7 @@ int vcs_Cantera_to_vprob(MultiPhase* mphase, VCS_PROB* vprob)
                 vprob->mf[kT] = mphase->moleFraction(kT);
              } else {
                throw CanteraError(" vcs_Cantera_to_vprob() ERROR",
-                                  "Unknown species type: " + int2str(vprob->SpeciesUnknownType[kT]));
+                                  "Unknown species type: {}", vprob->SpeciesUnknownType[kT]);
              }
 
             /*
diff --git a/src/equil/vcs_VolPhase.cpp b/src/equil/vcs_VolPhase.cpp
index b7c490fdd..8ff374465 100644
--- a/src/equil/vcs_VolPhase.cpp
+++ b/src/equil/vcs_VolPhase.cpp
@@ -397,8 +397,8 @@ void vcs_VolPhase::setMoleFractionsState(const double totalMoles,
     if (m_totalMolesInert > 0.0) {
         if (m_totalMolesInert > v_totalMoles) {
             throw CanteraError("vcs_VolPhase::setMolesFractionsState",
-                 "inerts greater than total: " + fp2str(v_totalMoles) + " " +
-                 fp2str(m_totalMolesInert));
+                 "inerts greater than total: {} {}",
+                 v_totalMoles, m_totalMolesInert);
         }
         fractotal = 1.0 - m_totalMolesInert/v_totalMoles;
     }
@@ -518,8 +518,7 @@ void vcs_VolPhase::setMolesFromVCSCheck(const int vcsStateStatus,
             Tcheck = v_totalMoles;
         } else {
             throw CanteraError("vcs_VolPhase::setMolesFromVCSCheck",
-                  "We have a consistency problem: " + fp2str(Tcheck) + " " +
-                  fp2str(v_totalMoles));
+                  "We have a consistency problem: {} {}", Tcheck, v_totalMoles);
         }
     }
 }
diff --git a/src/equil/vcs_nondim.cpp b/src/equil/vcs_nondim.cpp
index 606a1ec4a..894e0d8e1 100644
--- a/src/equil/vcs_nondim.cpp
+++ b/src/equil/vcs_nondim.cpp
@@ -30,8 +30,7 @@ double VCS_SOLVE::vcs_nondim_Farad(int mu_units, double TKelvin) const
     case VCS_UNITS_KELVIN:
         return ElectronCharge * Avogadro/ TKelvin;
     default:
-        throw CanteraError("vcs_nondim_Farad",
-                           "unknown units: " + int2str(mu_units));
+        throw CanteraError("vcs_nondim_Farad", "unknown units: {}", mu_units);
     }
 }
 
@@ -52,8 +51,7 @@ double VCS_SOLVE::vcs_nondimMult_TP(int mu_units, double TKelvin) const
     case VCS_UNITS_MKS:
         return TKelvin * GasConstant;
     default:
-        throw CanteraError("vcs_nondimMult_TP",
-                           "unknown units: " + int2str(mu_units));
+        throw CanteraError("vcs_nondimMult_TP", "unknown units: {}", mu_units);
     }
 }
 
@@ -101,8 +99,8 @@ void VCS_SOLVE::vcs_nondim_TP()
          */
         if (tmole_orig < 1.0E-200 || tmole_orig > 1.0E200) {
             throw CanteraError("VCS_SOLVE::vcs_nondim_TP",
-                               "Total input moles ," + fp2str(tmole_orig) +
-                               "is outside the range handled by vcs.\n");
+                "Total input moles, {} is outside the range handled by vcs.\n",
+                tmole_orig);
         }
 
         // Determine the scale of the problem
diff --git a/src/equil/vcs_solve.cpp b/src/equil/vcs_solve.cpp
index 8e6cb5335..de635dd2a 100644
--- a/src/equil/vcs_solve.cpp
+++ b/src/equil/vcs_solve.cpp
@@ -638,9 +638,9 @@ int VCS_SOLVE::vcs_prob_specifyFully(const VCS_PROB* pub)
             if (m_elemAbundancesGoal[i] != 0.0) {
                 if (fabs(m_elemAbundancesGoal[i]) > 1.0E-9) {
                     throw CanteraError("VCS_SOLVE::vcs_prob_specifyFully",
-                            "Charge neutrality condition " + m_elementName[i] +
-                            " is signicantly nonzero, " + fp2str(m_elemAbundancesGoal[i]) +
-                            ". Giving up");
+                            "Charge neutrality condition {} is signicantly "
+                            "nonzero, {}. Giving up",
+                            m_elementName[i], m_elemAbundancesGoal[i]);
                 } else {
                     if (m_debug_print_lvl >= 2) {
                         plogf("Charge neutrality condition %s not zero, %g. Setting it zero\n",
@@ -896,16 +896,15 @@ int VCS_SOLVE::vcs_prob_update(VCS_PROB* pub)
                 double tmp = m_molNumSpecies_old[k1];
                 if (! vcs_doubleEqual(pubPhase->electricPotential() , tmp)) {
                     throw CanteraError("VCS_SOLVE::vcs_prob_update",
-                            "We have an inconsistency in voltage, " +
-                            fp2str(pubPhase->electricPotential()) + " " +
-                            fp2str(tmp));
+                            "We have an inconsistency in voltage, {} {}",
+                            pubPhase->electricPotential(), tmp);
                 }
             }
 
             if (! vcs_doubleEqual(pub->mf[kT], vPhase->molefraction(k))) {
                 throw CanteraError("VCS_SOLVE::vcs_prob_update",
-                        "We have an inconsistency in mole fraction, " +
-                        fp2str(pub->mf[kT]) + " " + fp2str(vPhase->molefraction(k)));
+                        "We have an inconsistency in mole fraction, {} {}",
+                        pub->mf[kT], vPhase->molefraction(k));
             }
             if (pubPhase->speciesUnknownType(k) != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) {
                 sumMoles += pub->w[kT];
@@ -913,8 +912,8 @@ int VCS_SOLVE::vcs_prob_update(VCS_PROB* pub)
         }
         if (! vcs_doubleEqual(sumMoles, vPhase->totalMoles())) {
             throw CanteraError("VCS_SOLVE::vcs_prob_update",
-                            "We have an inconsistency in total moles, " +
-                            fp2str(sumMoles) + " " + fp2str(pubPhase->totalMoles()));
+                            "We have an inconsistency in total moles, {} {}",
+                            sumMoles, pubPhase->totalMoles());
         }
     }
 
diff --git a/src/equil/vcs_solve_TP.cpp b/src/equil/vcs_solve_TP.cpp
index c93d8c7db..30f501dc2 100644
--- a/src/equil/vcs_solve_TP.cpp
+++ b/src/equil/vcs_solve_TP.cpp
@@ -1066,9 +1066,8 @@ void VCS_SOLVE::solve_tp_inner(size_t& iti, size_t& it1,
             if (m_molNumSpecies_new[kspec] < 0.0 && (m_speciesUnknownType[kspec]
                     != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE)) {
                 throw CanteraError("VCS_SOLVE::solve_tp_inner",
-                    "vcs_solve_TP: ERROR on step change wt[" + int2str(kspec) + ":" +
-                    m_speciesName[kspec] + "]: " +
-                    fp2str(m_molNumSpecies_new[kspec]) + " < 0.0");
+                    "vcs_solve_TP: ERROR on step change wt[{}:{}]: {} < 0.0",
+                    kspec, m_speciesName[kspec], m_molNumSpecies_new[kspec]);
             }
         }
 
@@ -3272,8 +3271,7 @@ void VCS_SOLVE::vcs_dfe(const int stateCalc,
         molNum = &m_molNumSpecies_new[0];
     } else if (DEBUG_MODE_ENABLED) {
         throw CanteraError("VCS_SOLVE::vcs_dfe",
-                "Subroutine vcs_dfe called with bad stateCalc value: "+
-                int2str(stateCalc));
+                "Subroutine vcs_dfe called with bad stateCalc value: {}", stateCalc);
     }
 
     AssertThrowMsg(m_unitsState != VCS_DIMENSIONAL_G, "VCS_SOLVE::vcs_dfe",
@@ -3682,7 +3680,7 @@ void VCS_SOLVE::vcs_updateVP(const int vcsState)
                                          &m_tPhaseMoles_new[0]);
         } else if (DEBUG_MODE_ENABLED) {
             throw CanteraError("VCS_SOLVE::vcs_updateVP",
-                               "wrong stateCalc value: " + int2str(vcsState));
+                               "wrong stateCalc value: {}", vcsState);
         }
     }
 }
@@ -3742,7 +3740,7 @@ bool VCS_SOLVE::vcs_evaluate_speciesType()
                     break;
                 default:
                     throw CanteraError("VCS_SOLVE::vcs_evaluate_speciesType",
-                                       "Unknown type: " + int2str(m_speciesStatus[kspec]));
+                                       "Unknown type: {}", m_speciesStatus[kspec]);
                 }
             }
         }
diff --git a/src/equil/vcs_util.cpp b/src/equil/vcs_util.cpp
index 213f94322..99bdb9e64 100644
--- a/src/equil/vcs_util.cpp
+++ b/src/equil/vcs_util.cpp
@@ -73,8 +73,7 @@ double vcsUtil_gasConstant(int mu_units)
         /* joules / kg-mol K = kg m2 / s2 kg-mol K */
         return GasConstant;
     default:
-        throw CanteraError("vcsUtil_gasConstant",
-                           "uknown units: " + int2str(mu_units));
+        throw CanteraError("vcsUtil_gasConstant", "uknown units: {}", mu_units);
     }
 }
 
diff --git a/src/kinetics/AqueousKinetics.cpp b/src/kinetics/AqueousKinetics.cpp
index ba90f9792..f3617d96b 100644
--- a/src/kinetics/AqueousKinetics.cpp
+++ b/src/kinetics/AqueousKinetics.cpp
@@ -154,7 +154,7 @@ bool AqueousKinetics::addReaction(shared_ptr<Reaction> r)
         addElementaryReaction(dynamic_cast<ElementaryReaction&>(*r));
     } else {
         throw CanteraError("AqueousKinetics::addReaction",
-            "Invalid reaction type: " + int2str(r->reaction_type));
+            "Invalid reaction type: {}", r->reaction_type);
     }
     return true;
 }
diff --git a/src/kinetics/Falloff.cpp b/src/kinetics/Falloff.cpp
index d57769311..b24621ff2 100644
--- a/src/kinetics/Falloff.cpp
+++ b/src/kinetics/Falloff.cpp
@@ -14,8 +14,8 @@ void Falloff::init(const vector_fp& c)
 {
     if (c.size() != 0) {
         throw CanteraError("Falloff::init",
-            "Incorrect number of parameters. 0 required. Received " +
-            int2str(c.size()) + ".");
+            "Incorrect number of parameters. 0 required. Received {}.",
+            c.size());
     }
 }
 
@@ -23,8 +23,8 @@ void Troe::init(const vector_fp& c)
 {
     if (c.size() != 3 && c.size() != 4) {
         throw CanteraError("Troe::init",
-            "Incorrect number of parameters. 3 or 4 required. Received " +
-            int2str(c.size()) + ".");
+            "Incorrect number of parameters. 3 or 4 required. Received {}.",
+            c.size());
     }
     m_a = c[0];
     if (c[1] == 0.0) {
@@ -73,13 +73,13 @@ void SRI::init(const vector_fp& c)
 {
     if (c.size() != 3 && c.size() != 5) {
         throw CanteraError("SRI::init",
-            "Incorrect number of parameters. 3 or 5 required. Received " +
-            int2str(c.size()) + ".");
+            "Incorrect number of parameters. 3 or 5 required. Received {}.",
+            c.size());
     }
 
     if (c[2] < 0.0) {
         throw CanteraError("SRI::init()",
-                           "m_c parameter is less than zero: " + fp2str(c[2]));
+                           "m_c parameter is less than zero: {}", c[2]);
     }
     m_a = c[0];
     m_b = c[1];
@@ -88,7 +88,7 @@ void SRI::init(const vector_fp& c)
     if (c.size() == 5) {
         if (c[3] < 0.0) {
             throw CanteraError("SRI::init()",
-                               "m_d parameter is less than zero: " + fp2str(c[3]));
+                               "m_d parameter is less than zero: {}", c[3]);
         }
         m_d = c[3];
         m_e = c[4];
diff --git a/src/kinetics/GasKinetics.cpp b/src/kinetics/GasKinetics.cpp
index c409818ef..8ac2a0a70 100644
--- a/src/kinetics/GasKinetics.cpp
+++ b/src/kinetics/GasKinetics.cpp
@@ -264,7 +264,7 @@ bool GasKinetics::addReaction(shared_ptr<Reaction> r)
         break;
     default:
         throw CanteraError("GasKinetics::addReaction",
-            "Unknown reaction type specified: " + int2str(r->reaction_type));
+            "Unknown reaction type specified: {}", r->reaction_type);
     }
     return true;
 }
@@ -356,7 +356,7 @@ void GasKinetics::modifyReaction(size_t i, shared_ptr<Reaction> rNew)
         break;
     default:
         throw CanteraError("GasKinetics::modifyReaction",
-            "Unknown reaction type specified: " + int2str(rNew->reaction_type));
+            "Unknown reaction type specified: {}", rNew->reaction_type);
     }
 
     // invalidate all cached data
diff --git a/src/kinetics/InterfaceKinetics.cpp b/src/kinetics/InterfaceKinetics.cpp
index 0d339b3e0..f88a7381a 100644
--- a/src/kinetics/InterfaceKinetics.cpp
+++ b/src/kinetics/InterfaceKinetics.cpp
@@ -218,7 +218,7 @@ void InterfaceKinetics::updateKc()
         for (size_t i = 0; i < m_nrev; i++) {
             size_t irxn = m_revindex[i];
             if (irxn == npos || irxn >= nReactions()) {
-                throw CanteraError("InterfaceKinetics", "illegal value: irxn = "+int2str(irxn));
+                throw CanteraError("InterfaceKinetics", "illegal value: irxn = {}", irxn);
             }
             // WARNING this may overflow HKM
             m_rkcn[irxn] = exp(m_rkcn[irxn]*rrt);
@@ -892,8 +892,8 @@ void InterfaceKinetics::finalize()
     m_surf = (SurfPhase*)&thermo(ks);
     if (m_surf->nDim() != 2) {
         throw CanteraError("InterfaceKinetics::finalize",
-                           "expected interface dimension = 2, but got dimension = "
-                           +int2str(m_surf->nDim()));
+                           "expected interface dimension = 2, but got dimension = {}",
+                           m_surf->nDim());
     }
     m_StandardConc.resize(m_kk, 0.0);
     m_deltaG0.resize(nReactions(), 0.0);
@@ -1077,8 +1077,8 @@ void EdgeKinetics::finalize()
     m_surf = (SurfPhase*)&thermo(ks);
     if (m_surf->nDim() != 1) {
         throw CanteraError("EdgeKinetics::finalize",
-                           "expected interface dimension = 1, but got dimension = "
-                           +int2str(m_surf->nDim()));
+                           "expected interface dimension = 1, but got dimension = {}",
+                           m_surf->nDim());
     }
     m_StandardConc.resize(m_kk, 0.0);
     m_deltaG0.resize(safe_reaction_size, 0.0);
diff --git a/src/kinetics/Kinetics.cpp b/src/kinetics/Kinetics.cpp
index 72179eda2..8735ec83e 100644
--- a/src/kinetics/Kinetics.cpp
+++ b/src/kinetics/Kinetics.cpp
@@ -221,10 +221,10 @@ std::pair<size_t, size_t> Kinetics::checkDuplicates(bool throw_err) const
                 }
             }
             if (throw_err) {
-                string msg = string("Undeclared duplicate reactions detected:\n")
-                             +"Reaction "+int2str(i+1)+": "+other.equation()
-                             +"\nReaction "+int2str(m+1)+": "+R.equation()+"\n";
-                throw CanteraError("installReaction", msg);
+                throw CanteraError("installReaction",
+                        "Undeclared duplicate reactions detected:\n"
+                        "Reaction {}: {}\nReaction {}: {}\n",
+                        int2str(i+1), other.equation(), int2str(m+1), R.equation());
             } else {
                 return make_pair(i,m);
             }
@@ -304,8 +304,8 @@ void Kinetics::checkReactionBalance(const Reaction& R)
         double elemdiff = fabs(balp[elem] - balr[elem]);
         if (elemsum > 0.0 && elemdiff/elemsum > 1e-4) {
             ok = false;
-            msg += "  " + elem + "           " + fp2str(balr[elem]) +
-                   "           " + fp2str(balp[elem]) + "\n";
+            msg += fmt::format("  {}           {}           {}\n",
+                               elem, balr[elem], balp[elem]);
         }
     }
     if (!ok) {
@@ -394,7 +394,7 @@ size_t Kinetics::speciesPhaseIndex(size_t k)
             return n;
         }
     }
-    throw CanteraError("speciesPhaseIndex", "illegal species index: "+int2str(k));
+    throw CanteraError("speciesPhaseIndex", "illegal species index: {}", k);
     return npos;
 }
 
@@ -632,20 +632,20 @@ void Kinetics::modifyReaction(size_t i, shared_ptr<Reaction> rNew)
     shared_ptr<Reaction>& rOld = m_reactions[i];
     if (rNew->reaction_type != rOld->reaction_type) {
         throw CanteraError("Kinetics::modifyReaction",
-            "Reaction types are different: " + int2str(rOld->reaction_type) +
-            " != " + int2str(rNew->reaction_type) + ".");
+            "Reaction types are different: {} != {}.",
+            rOld->reaction_type, rNew->reaction_type);
     }
 
     if (rNew->reactants != rOld->reactants) {
         throw CanteraError("Kinetics::modifyReaction",
-            "Reactants are different: '" + rOld->reactantString() + "' != '" +
-            rNew->reactantString() + "'.");
+            "Reactants are different: '{}' != '{}'.",
+            rOld->reactantString(), rNew->reactantString());
     }
 
     if (rNew->products != rOld->products) {
         throw CanteraError("Kinetics::modifyReaction",
-            "Products are different: '" + rOld->productString() + "' != '" +
-            rNew->productString() + "'.");
+            "Products are different: '{}' != '{}'.",
+            rOld->productString(), rNew->productString());
     }
     m_reactions[i] = rNew;
 }
diff --git a/src/kinetics/Reaction.cpp b/src/kinetics/Reaction.cpp
index f8fd79988..4f28506f7 100644
--- a/src/kinetics/Reaction.cpp
+++ b/src/kinetics/Reaction.cpp
@@ -293,23 +293,23 @@ void readFalloff(FalloffReaction& R, const XML_Node& rc_node)
         falloff_type = SIMPLE_FALLOFF;
         if (np != 0) {
             throw CanteraError("readFalloff", "Lindemann parameterization "
-                "takes no parameters, but " + int2str(np) + "were given");
+                "takes no parameters, but {} were given", np);
         }
     } else if (lowercase(falloff["type"]) == "troe") {
         falloff_type = TROE_FALLOFF;
         if (np != 3 && np != 4) {
             throw CanteraError("readFalloff", "Troe parameterization takes "
-                "3 or 4 parameters, but " + int2str(np) + "were given");
+                "3 or 4 parameters, but {} were given", np);
         }
     } else if (lowercase(falloff["type"]) == "sri") {
         falloff_type = SRI_FALLOFF;
         if (np != 3 && np != 5) {
             throw CanteraError("readFalloff", "SRI parameterization takes "
-                "3 or 5 parameters, but " + int2str(np) + "were given");
+                "3 or 5 parameters, but {} were given", np);
         }
     } else {
-        throw CanteraError("readFalloff", "Unrecognized falloff type: '" +
-            falloff["type"] + "'");
+        throw CanteraError("readFalloff", "Unrecognized falloff type: '{}'",
+                           falloff["type"]);
     }
     R.falloff = newFalloff(falloff_type, falloff_parameters);
 }
diff --git a/src/kinetics/RxnRates.cpp b/src/kinetics/RxnRates.cpp
index a04e7ab23..c4ad7416e 100644
--- a/src/kinetics/RxnRates.cpp
+++ b/src/kinetics/RxnRates.cpp
@@ -105,9 +105,8 @@ void Plog::validate(const std::string& equation)
                 // message will correctly indicate that the problematic rate
                 // expression is at the higher of the adjacent pressures.
                 throw CanteraError("Plog::validate",
-                        "Invalid rate coefficient for reaction '" + equation +
-                        "'\nat P = " + fp2str(std::exp((++iter)->first)) +
-                        ", T = " + fp2str(T[i]));
+                    "Invalid rate coefficient for reaction '{}'\nat P = {}, T = {}",
+                    equation, std::exp((++iter)->first), T[i]);
             }
         }
     }
diff --git a/src/kinetics/solveSP.cpp b/src/kinetics/solveSP.cpp
index 3ff744883..75685b7d9 100644
--- a/src/kinetics/solveSP.cpp
+++ b/src/kinetics/solveSP.cpp
@@ -787,7 +787,7 @@ void solveSP::print_header(int ioflag, int ifunc, doublereal time_scale,
             writelogf("           for a total of %9.3e sec\n", time_scale);
         } else {
             throw CanteraError("solveSP::print_header",
-                               "Unknown ifunc flag = " + int2str(ifunc));
+                               "Unknown ifunc flag = {}", ifunc);
         }
 
         if (m_bulkFunc == BULK_DEPOSITION) {
@@ -796,7 +796,7 @@ void solveSP::print_header(int ioflag, int ifunc, doublereal time_scale,
             writelog("     Bulk Phases have fixed compositions\n");
         } else {
             throw CanteraError("solveSP::print_header",
-                               "Unknown bulkFunc flag = " + int2str(m_bulkFunc));
+                               "Unknown bulkFunc flag = {}", m_bulkFunc);
         }
 
         if (damping) {
diff --git a/src/numerics/BandMatrix.cpp b/src/numerics/BandMatrix.cpp
index 9199cba05..b9684cbcb 100644
--- a/src/numerics/BandMatrix.cpp
+++ b/src/numerics/BandMatrix.cpp
@@ -338,7 +338,7 @@ doublereal BandMatrix::rcond(doublereal a1norm)
             writelogf("BandMatrix::rcond(): DGBCON returned INFO = %d\n", rinfo);
         }
         if (! useReturnErrorCode) {
-            throw CanteraError("BandMatrix::rcond()", "DGBCON returned INFO = " + int2str(rinfo));
+            throw CanteraError("BandMatrix::rcond()", "DGBCON returned INFO = {}", rinfo);
         }
     }
     return rcond;
diff --git a/src/numerics/polyfit.cpp b/src/numerics/polyfit.cpp
index b3d53390e..9435a8298 100644
--- a/src/numerics/polyfit.cpp
+++ b/src/numerics/polyfit.cpp
@@ -40,9 +40,8 @@ doublereal polyfit(int n, doublereal* x, doublereal* y, doublereal* w,
              &ierr, &awork[0]);
     if (ierr != 1) {
         throw CanteraError("polyfit",
-                           "DPOLFT returned error code IERR = " + int2str(ierr) +
-                           "while attempting to fit " + int2str(n) + " data points "
-                           + "to a polynomial of degree " + int2str(maxdeg));
+            "DPOLFT returned error code IERR = {} while attempting to fit {}"
+            " data points  to a polynomial of degree {}", ierr, n, maxdeg);
     }
     ndeg = ndg;
     _DPCOEF_(&ndg, &zer, r, &awork[0]);
diff --git a/src/oneD/MultiNewton.cpp b/src/oneD/MultiNewton.cpp
index f470cf451..091717d06 100644
--- a/src/oneD/MultiNewton.cpp
+++ b/src/oneD/MultiNewton.cpp
@@ -196,14 +196,11 @@ void MultiNewton::step(doublereal* x, doublereal* step,
         size_t pt = offset/dom.nComponents();
         size_t comp = offset - pt*dom.nComponents();
         throw CanteraError("MultiNewton::step",
-                           "Jacobian is singular for domain "+
-                           dom.id() + ", component "
-                           +dom.componentName(comp)+" at point "
-                           +int2str(pt)+"\n(Matrix row "
-                           +int2str(iok)+") \nsee file bandmatrix.csv\n");
+            "Jacobian is singular for domain {}, component {} at point {}\n"
+            "(Matrix row {}) \nsee file bandmatrix.csv\n",
+            dom.id(), dom.componentName(comp), pt, iok);
     } else if (int(iok) < 0) {
-        throw CanteraError("MultiNewton::step",
-                           "iok = "+int2str(iok));
+        throw CanteraError("MultiNewton::step", "iok = {}", iok);
     }
 }
 
diff --git a/src/oneD/Sim1D.cpp b/src/oneD/Sim1D.cpp
index a531767ed..1fc036174 100644
--- a/src/oneD/Sim1D.cpp
+++ b/src/oneD/Sim1D.cpp
@@ -124,9 +124,8 @@ void Sim1D::restore(const std::string& fname, const std::string& id,
     vector<XML_Node*> xd = f->getChildren("domain");
     if (xd.size() != m_nd) {
         throw CanteraError("Sim1D::restore", "Solution does not contain the "
-                           " correct number of domains. Found " +
-                           int2str(xd.size()) + "expected " +
-                           int2str(m_nd) + ".\n");
+            " correct number of domains. Found {} expected {}.\n",
+            xd.size(), m_nd);
     }
     size_t sz = 0;
     for (size_t m = 0; m < m_nd; m++) {
@@ -207,7 +206,7 @@ int Sim1D::newtonSolve(int loglevel)
         return -1;
     } else {
         throw CanteraError("Sim1D::newtonSolve",
-                           "ERROR: OneDim::solve returned m = " + int2str(m) + "\n");
+                           "ERROR: OneDim::solve returned m = {}", m);
     }
 }
 
@@ -236,7 +235,7 @@ void Sim1D::solve(int loglevel, bool refine_grid)
                     writelog("    success.\n\n");
                     writelog("Problem solved on [");
                     for (size_t mm = 1; mm < nDomains(); mm+=2) {
-                        writelog(int2str(domain(mm).nPoints()));
+                        writelog("{}", domain(mm).nPoints());
                         if (mm + 2 < nDomains()) {
                             writelog(", ");
                         }
diff --git a/src/oneD/StFlow.cpp b/src/oneD/StFlow.cpp
index e2bd00113..065e1adbd 100644
--- a/src/oneD/StFlow.cpp
+++ b/src/oneD/StFlow.cpp
@@ -754,9 +754,8 @@ void StFlow::restore(const XML_Node& dom, doublereal* soln, int loglevel)
                 m_do_energy[i] = x[i];
             }
         } else if (!x.empty()) {
-            throw CanteraError("StFlow::restore", "energy_enabled is length" +
-                               int2str(x.size()) + "but should be length" +
-                               int2str(nPoints()));
+            throw CanteraError("StFlow::restore", "energy_enabled is length {}"
+                               "but should be length {}", x.size(), nPoints());
         }
     }
 
@@ -770,9 +769,9 @@ void StFlow::restore(const XML_Node& dom, doublereal* soln, int loglevel)
             // This may occur when restoring from a mechanism with a different
             // number of species.
             if (loglevel > 0) {
-                writelog("\nWarning: StFlow::restore: species_enabled is length " +
-                         int2str(x.size()) + " but should be length " +
-                         int2str(m_nsp) + ". Enabling all species equations by default.");
+                writelog("\nWarning: StFlow::restore: species_enabled is "
+                    "length {} but should be length {}. Enabling all species "
+                    "equations by default.", x.size(), m_nsp);
             }
             m_do_species.assign(m_nsp, true);
         }
diff --git a/src/oneD/boundaries1D.cpp b/src/oneD/boundaries1D.cpp
index 8d3d1c96d..3097a3985 100644
--- a/src/oneD/boundaries1D.cpp
+++ b/src/oneD/boundaries1D.cpp
@@ -50,9 +50,8 @@ void Bdry1D::_init(size_t n)
             m_phase_left = &m_flow_left->phase();
         } else {
             throw CanteraError("Bdry1D::_init",
-                               "Boundary domains can only be "
-                               "connected on the left to flow domains, not type "+int2str(r.domainType())
-                               + " domains.");
+                "Boundary domains can only be connected on the left to flow "
+                "domains, not type {} domains.", r.domainType());
         }
     }
 
@@ -68,9 +67,8 @@ void Bdry1D::_init(size_t n)
             m_phase_right = &m_flow_right->phase();
         } else {
             throw CanteraError("Bdry1D::_init",
-                               "Boundary domains can only be "
-                               "connected on the right to flow domains, not type "+int2str(r.domainType())
-                               + " domains.");
+                "Boundary domains can only be connected on the right to flow "
+                "domains, not type {} domains.", r.domainType());
         }
     }
 }
diff --git a/src/oneD/refine.cpp b/src/oneD/refine.cpp
index cad46fee9..fd040f3ad 100644
--- a/src/oneD/refine.cpp
+++ b/src/oneD/refine.cpp
@@ -21,20 +21,17 @@ void Refiner::setCriteria(doublereal ratio, doublereal slope,
 {
     if (ratio < 2.0) {
         throw CanteraError("Refiner::setCriteria",
-            "'ratio' must be greater than 2.0 (" + fp2str(ratio) +
-            " was specified).");
+            "'ratio' must be greater than 2.0 ({} was specified).", ratio);
     } else if (slope < 0.0 || slope > 1.0) {
         throw CanteraError("Refiner::setCriteria",
-            "'slope' must be between 0.0 and 1.0 (" + fp2str(slope) +
-            " was specified).");
+            "'slope' must be between 0.0 and 1.0 ({} was specified).", slope);
     } else if (curve < 0.0 || curve > 1.0) {
         throw CanteraError("Refiner::setCriteria",
-            "'curve' must be between 0.0 and 1.0 (" + fp2str(curve) +
-            " was specified).");
+            "'curve' must be between 0.0 and 1.0 ({} was specified).", curve);
     } else if (prune > curve || prune > slope) {
         throw CanteraError("Refiner::setCriteria",
-            "'prune' must be less than 'curve' and 'slope' (" + fp2str(prune) +
-            " was specified).");
+            "'prune' must be less than 'curve' and 'slope' ({} was specified).",
+            prune);
     }
     m_ratio = ratio;
     m_slope = slope;
@@ -46,7 +43,7 @@ int Refiner::analyze(size_t n, const doublereal* z,
                      const doublereal* x)
 {
     if (n >= m_npmax) {
-        writelog("max number of grid points reached ("+int2str(m_npmax)+".\n");
+        writelog("max number of grid points reached ({}).\n", m_npmax);
         return -2;
     }
 
@@ -214,7 +211,7 @@ void Refiner::show()
                  m_domain->id()+".\n"
                  +"    New points inserted after grid points ");
         for (const auto& loc : m_loc) {
-            writelog(int2str(loc.first)+" ");
+            writelog("{} ", loc.first);
         }
         writelog("\n");
         writelog("    to resolve ");
diff --git a/src/thermo/GibbsExcessVPSSTP.cpp b/src/thermo/GibbsExcessVPSSTP.cpp
index 18484f9ea..15aa05fc6 100644
--- a/src/thermo/GibbsExcessVPSSTP.cpp
+++ b/src/thermo/GibbsExcessVPSSTP.cpp
@@ -192,7 +192,7 @@ double GibbsExcessVPSSTP::checkMFSum(const doublereal* const x) const
     doublereal norm = std::accumulate(x, x + m_kk, 0.0);
     if (fabs(norm - 1.0) > 1.0E-9) {
         throw CanteraError("GibbsExcessVPSSTP::checkMFSum",
-                           "(MF sum - 1) exceeded tolerance of 1.0E-9:" + fp2str(norm));
+            "(MF sum - 1) exceeded tolerance of 1.0E-9: {}", norm);
     }
     return norm;
 }
diff --git a/src/thermo/IonsFromNeutralVPSSTP.cpp b/src/thermo/IonsFromNeutralVPSSTP.cpp
index 4e10bbd96..4654177e0 100644
--- a/src/thermo/IonsFromNeutralVPSSTP.cpp
+++ b/src/thermo/IonsFromNeutralVPSSTP.cpp
@@ -519,7 +519,7 @@ void IonsFromNeutralVPSSTP::calcNeutralMoleculeMoleFractions() const
         }
         if (fabs(sum) > 1.0E-11) {
             throw CanteraError("IonsFromNeutralVPSSTP::calcNeutralMoleculeMoleFractions",
-                               "molefracts don't sum to one: " + fp2str(sum));
+                               "molefracts don't sum to one: {}", sum);
         }
     }
 
diff --git a/src/thermo/MargulesVPSSTP.cpp b/src/thermo/MargulesVPSSTP.cpp
index b10b6f6ed..95f6d3ab1 100644
--- a/src/thermo/MargulesVPSSTP.cpp
+++ b/src/thermo/MargulesVPSSTP.cpp
@@ -608,7 +608,8 @@ void MargulesVPSSTP::readXMLBinarySpecies(XML_Node& xmLBinarySpecies)
     //   @TODO Figure out what the original reason is for putting an error condition for charged species
     //         Seems OK to me.
     if (charge(aSpecies) != 0.0) {
-        throw CanteraError("MargulesVPSSTP::readXMLBinarySpecies", "speciesA has a charge: " + fp2str(charge(aSpecies)));
+        throw CanteraError("MargulesVPSSTP::readXMLBinarySpecies",
+            "speciesA has a charge: {}", charge(aSpecies));
     }
     size_t bSpecies = speciesIndex(bName);
     if (bSpecies == npos) {
@@ -616,7 +617,8 @@ void MargulesVPSSTP::readXMLBinarySpecies(XML_Node& xmLBinarySpecies)
     }
     string bspName = speciesName(bSpecies);
     if (charge(bSpecies) != 0.0) {
-        throw CanteraError("MargulesVPSSTP::readXMLBinarySpecies", "speciesB has a charge: " + fp2str(charge(bSpecies)));
+        throw CanteraError("MargulesVPSSTP::readXMLBinarySpecies",
+            "speciesB has a charge: {}", charge(bSpecies));
     }
 
     resizeNumInteractions(numBinaryInteractions_ + 1);
diff --git a/src/thermo/MineralEQ3.cpp b/src/thermo/MineralEQ3.cpp
index 1594d596f..74998f272 100644
--- a/src/thermo/MineralEQ3.cpp
+++ b/src/thermo/MineralEQ3.cpp
@@ -283,8 +283,8 @@ void MineralEQ3::convertDGFormation()
     // If the discrepancy is greater than 100 cal gmol-1, print an error
     if (fabs(Hcalc -DHjmol) > 10.* 1.0E6 * 4.184) {
         throw CanteraError("installMinEQ3asShomateThermoFromXML()",
-                           "DHjmol is not consistent with G and S" +
-                           fp2str(Hcalc) + " vs " + fp2str(DHjmol));
+                           "DHjmol is not consistent with G and S: {} vs {}",
+                           Hcalc, DHjmol);
     }
 }
 
diff --git a/src/thermo/MixtureFugacityTP.cpp b/src/thermo/MixtureFugacityTP.cpp
index a2ad047fa..ebc71450e 100644
--- a/src/thermo/MixtureFugacityTP.cpp
+++ b/src/thermo/MixtureFugacityTP.cpp
@@ -660,8 +660,8 @@ int MixtureFugacityTP::corr0(doublereal TKelvin, doublereal pres, doublereal& de
     if (densGas <= 0.0) {
         if (retn == -1) {
             throw CanteraError("MixtureFugacityTP::corr0",
-                               "Error occurred trying to find gas density at (T,P) = "
-                               + fp2str(TKelvin) + "  " + fp2str(pres));
+                "Error occurred trying to find gas density at (T,P) = {}  {}",
+                TKelvin, pres);
         }
         retn = -2;
     } else {
diff --git a/src/thermo/MolalityVPSSTP.cpp b/src/thermo/MolalityVPSSTP.cpp
index c43a5ccb0..93d442cd9 100644
--- a/src/thermo/MolalityVPSSTP.cpp
+++ b/src/thermo/MolalityVPSSTP.cpp
@@ -84,7 +84,7 @@ void MolalityVPSSTP::setpHScale(const int pHscaleType)
     m_pHScalingType = pHscaleType;
     if (pHscaleType != PHSCALE_PITZER && pHscaleType != PHSCALE_NBS) {
         throw CanteraError("MolalityVPSSTP::setpHScale",
-                           "Unknown scale type: " + int2str(pHscaleType));
+                           "Unknown scale type: {}", pHscaleType);
     }
 }
 
diff --git a/src/thermo/NasaPoly2.cpp b/src/thermo/NasaPoly2.cpp
index e34aa1642..8222032b1 100644
--- a/src/thermo/NasaPoly2.cpp
+++ b/src/thermo/NasaPoly2.cpp
@@ -13,37 +13,28 @@ void NasaPoly2::validate(const std::string& name)
 
     double delta = cp_low - cp_high;
     if (fabs(delta/(fabs(cp_low)+1.0E-4)) > 0.001) {
-        writelog("\n\n**** WARNING ****\nFor species "+name+
-                 ", discontinuity in cp/R detected at Tmid = "
-                 +fp2str(m_midT)+"\n");
-        writelog("\tValue computed using low-temperature polynomial:  "
-                 +fp2str(cp_low)+".\n");
-        writelog("\tValue computed using high-temperature polynomial: "
-                 +fp2str(cp_high)+".\n");
+        writelog("\n\n**** WARNING ****\nFor species {}, discontinuity"
+                 " in cp/R detected at Tmid = {}\n", name, m_midT);
+        writelog("\tValue computed using low-temperature polynomial:  {}\n", cp_low);
+        writelog("\tValue computed using high-temperature polynomial: {}\n", cp_high);
     }
 
     // enthalpy
     delta = h_low - h_high;
     if (fabs(delta/(fabs(h_low)+cp_low*m_midT)) > 0.001) {
-        writelog("\n\n**** WARNING ****\nFor species "+name+
-                 ", discontinuity in h/RT detected at Tmid = "
-                 +fp2str(m_midT)+"\n");
-        writelog("\tValue computed using low-temperature polynomial:  "
-                 +fp2str(h_low)+".\n");
-        writelog("\tValue computed using high-temperature polynomial: "
-                 +fp2str(h_high)+".\n");
+        writelog("\n\n**** WARNING ****\nFor species {}, discontinuity"
+                 " in h/RT detected at Tmid = {}\n", name, m_midT);
+        writelog("\tValue computed using low-temperature polynomial:  {}\n", h_low);
+        writelog("\tValue computed using high-temperature polynomial: {}\n", h_high);
     }
 
     // entropy
     delta = s_low - s_high;
     if (fabs(delta/(fabs(s_low)+cp_low)) > 0.001) {
-        writelog("\n\n**** WARNING ****\nFor species "+name+
-                 ", discontinuity in s/R detected at Tmid = "
-                 +fp2str(m_midT)+"\n");
-        writelog("\tValue computed using low-temperature polynomial:  "
-                 +fp2str(s_low)+".\n");
-        writelog("\tValue computed using high-temperature polynomial: "
-                 +fp2str(s_high)+".\n");
+        writelog("\n\n**** WARNING ****\nFor species {}, discontinuity"
+                 " in s/R detected at Tmid = {}\n", name, m_midT);
+        writelog("\tValue computed using low-temperature polynomial:  {}\n", s_low);
+        writelog("\tValue computed using high-temperature polynomial: {}\n", s_high);
     }
 }
 
diff --git a/src/thermo/PDSS_HKFT.cpp b/src/thermo/PDSS_HKFT.cpp
index f09539811..6832d2466 100644
--- a/src/thermo/PDSS_HKFT.cpp
+++ b/src/thermo/PDSS_HKFT.cpp
@@ -430,16 +430,14 @@ void PDSS_HKFT::initThermo()
     if (fabs(Hcalc -DHjmol) > 100.* 1.0E3 * 4.184) {
         std::string sname = m_tp->speciesName(m_spindex);
         if (s_InputInconsistencyErrorExit) {
-            throw CanteraError(" PDSS_HKFT::initThermo() for " + sname,
-                               "DHjmol is not consistent with G and S: " +
-                               fp2str(Hcalc/(4.184E3)) + " vs "
-                               + fp2str(m_deltaH_formation_tr_pr) + "cal gmol-1");
+            throw CanteraError("PDSS_HKFT::initThermo()", "For {}, DHjmol is"
+                " not consistent with G and S: {} vs {} cal gmol-1",
+                sname, Hcalc/4.184E3, m_deltaH_formation_tr_pr);
         } else {
-            writelog(" PDSS_HKFT::initThermo() WARNING: "
-                     "DHjmol for " + sname + " is not consistent with G and S: calculated " +
-                     fp2str(Hcalc/(4.184E3)) + " vs input "
-                     + fp2str(m_deltaH_formation_tr_pr) + "cal gmol-1");
-            writelog("                                : continuing with consistent DHjmol = " + fp2str(Hcalc/(4.184E3)));
+            writelog("PDSS_HKFT::initThermo() WARNING: DHjmol for {} is not"
+                " consistent with G and S: calculated {} vs input {} cal gmol-1",
+                sname, Hcalc/4.184E3, m_deltaH_formation_tr_pr);
+            writelog("                                : continuing with consistent DHjmol = {}", Hcalc/4.184E3);
             m_deltaH_formation_tr_pr = Hcalc / (1.0E3 * 4.184);
         }
     }
diff --git a/src/thermo/PDSS_Water.cpp b/src/thermo/PDSS_Water.cpp
index f9124118c..a46a5efaf 100644
--- a/src/thermo/PDSS_Water.cpp
+++ b/src/thermo/PDSS_Water.cpp
@@ -310,10 +310,8 @@ void PDSS_Water::setPressure(doublereal p)
 
     doublereal dd = m_sub.density(T, p, waterState, dens);
     if (dd <= 0.0) {
-        std::string stateString = "T = " +
-                                  fp2str(T) + " K and p = " + fp2str(p) + " Pa";
         throw CanteraError("PDSS_Water:setPressure()",
-                           "Failed to set water SS state: " + stateString);
+            "Failed to set water SS state: T = {} K and p = {} Pa", T, p);
     }
     m_dens = dd;
     m_pres = p;
@@ -339,7 +337,7 @@ doublereal PDSS_Water::dthermalExpansionCoeffdT() const
     doublereal dd = m_sub.density(tt, pres, m_iState, m_dens);
     if (dd < 0.0) {
         throw CanteraError("PDSS_Water::dthermalExpansionCoeffdT",
-                           "unable to solve for the density at T = " + fp2str(tt) + ", P = " + fp2str(pres));
+            "unable to solve for the density at T = {}, P = {}", tt, pres);
     }
     doublereal vald = m_sub.coeffThermExp();
     m_sub.setState_TR(m_temp, dens_save);
diff --git a/src/thermo/RedlichKwongMFTP.cpp b/src/thermo/RedlichKwongMFTP.cpp
index 7a88e1174..f02845dfe 100644
--- a/src/thermo/RedlichKwongMFTP.cpp
+++ b/src/thermo/RedlichKwongMFTP.cpp
@@ -1289,10 +1289,9 @@ int RedlichKwongMFTP::NicholsSolve(double TKelvin, double pres, doublereal a, do
             if (fabs(tmp) > 1.0E-4) {
                 for (int j = 0; j < 3; j++) {
                     if (j != i && fabs(Vroot[i] - Vroot[j]) < 1.0E-4 * (fabs(Vroot[i]) + fabs(Vroot[j]))) {
-                        writelog("RedlichKwongMFTP::NicholsSolve(T = " + fp2str(TKelvin) + ", p = " +
-                                 fp2str(pres) + "): WARNING roots have merged: " +
-                                 fp2str(Vroot[i]) + ", " + fp2str(Vroot[j]));
-                        writelogendl();
+                        writelog("RedlichKwongMFTP::NicholsSolve(T = {}, p = {}):"
+                                 " WARNING roots have merged: {}, {}\n",
+                                 TKelvin, pres, Vroot[i], Vroot[j]);
                     }
                 }
             }
@@ -1351,8 +1350,8 @@ int RedlichKwongMFTP::NicholsSolve(double TKelvin, double pres, doublereal a, do
             }
         }
         if ((fabs(res) > 1.0E-14) && (fabs(res) > 1.0E-14 * fabs(dresdV) * fabs(Vroot[i]))) {
-            writelog("RedlichKwongMFTP::NicholsSolve(T = " + fp2str(TKelvin) + ", p = " +
-                     fp2str(pres) + "): WARNING root didn't converge V = " + fp2str(Vroot[i]));
+            writelog("RedlichKwongMFTP::NicholsSolve(T = {}, p = {}): "
+                "WARNING root didn't converge V = {}", TKelvin, pres, Vroot[i]);
             writelogendl();
         }
     }
diff --git a/src/thermo/SingleSpeciesTP.cpp b/src/thermo/SingleSpeciesTP.cpp
index 9782ea53b..fba5a30ce 100644
--- a/src/thermo/SingleSpeciesTP.cpp
+++ b/src/thermo/SingleSpeciesTP.cpp
@@ -240,7 +240,7 @@ void SingleSpeciesTP::setState_HP(doublereal h, doublereal p,
             return;
         }
     }
-    throw CanteraError("setState_HP","no convergence. dt = " + fp2str(dt));
+    throw CanteraError("setState_HP","no convergence. dt = {}", dt);
 }
 
 void SingleSpeciesTP::setState_UV(doublereal u, doublereal v,
@@ -259,9 +259,8 @@ void SingleSpeciesTP::setState_UV(doublereal u, doublereal v,
             return;
         }
     }
-    throw CanteraError("setState_UV",
-                       "no convergence. dt = " + fp2str(dt)+"\n"
-                       +"u = "+fp2str(u)+" v = "+fp2str(v)+"\n");
+    throw CanteraError("setState_UV", "no convergence. dt = {}\n"
+                       "u = {} v = {}\n", dt, u, v);
 }
 
 void SingleSpeciesTP::setState_SP(doublereal s, doublereal p,
@@ -276,7 +275,7 @@ void SingleSpeciesTP::setState_SP(doublereal s, doublereal p,
             return;
         }
     }
-    throw CanteraError("setState_SP","no convergence. dt = " + fp2str(dt));
+    throw CanteraError("setState_SP","no convergence. dt = {}", dt);
 }
 
 void SingleSpeciesTP::setState_SV(doublereal s, doublereal v,
@@ -295,7 +294,7 @@ void SingleSpeciesTP::setState_SV(doublereal s, doublereal v,
             return;
         }
     }
-    throw CanteraError("setState_SV","no convergence. dt = " + fp2str(dt));
+    throw CanteraError("setState_SV","no convergence. dt = {}", dt);
 }
 
 void SingleSpeciesTP::initThermo()
diff --git a/src/thermo/SpeciesThermoFactory.cpp b/src/thermo/SpeciesThermoFactory.cpp
index 53eac0431..b093082b3 100644
--- a/src/thermo/SpeciesThermoFactory.cpp
+++ b/src/thermo/SpeciesThermoFactory.cpp
@@ -48,7 +48,7 @@ SpeciesThermoInterpType* newSpeciesThermoInterpType(int type, double tlow,
         return new Adsorbate(tlow, thigh, pref, coeffs);
     default:
         throw CanteraError("newSpeciesThermoInterpType",
-                           "Unknown species thermo type: " + int2str(type) + ".");
+                           "Unknown species thermo type: {}.", type);
     }
 }
 
diff --git a/src/thermo/SurfPhase.cpp b/src/thermo/SurfPhase.cpp
index 2ad091b3e..e233ed898 100644
--- a/src/thermo/SurfPhase.cpp
+++ b/src/thermo/SurfPhase.cpp
@@ -254,7 +254,7 @@ void SurfPhase::setSiteDensity(doublereal n0)
 {
     if (n0 <= 0.0) {
         throw CanteraError("SurfPhase::setSiteDensity",
-                           "Site density must be positive. Got " + fp2str(n0));
+                           "Site density must be positive. Got {}", n0);
     }
     m_n0 = n0;
     m_logn0 = log(m_n0);
diff --git a/src/thermo/ThermoFactory.cpp b/src/thermo/ThermoFactory.cpp
index 26a428e56..24abb53c6 100644
--- a/src/thermo/ThermoFactory.cpp
+++ b/src/thermo/ThermoFactory.cpp
@@ -450,8 +450,8 @@ void importPhase(XML_Node& phase, ThermoPhase* th)
 
     size_t nsp = spDataNodeList.size();
     if (ssConvention == cSS_CONVENTION_SLAVE && nsp > 0) {
-        throw CanteraError("importPhase()", "For Slave standard states, number of species must be zero: "
-                           + int2str(nsp));
+        throw CanteraError("importPhase()", "For Slave standard states, "
+            "number of species must be zero: {}", nsp);
     }
     for (size_t k = 0; k < nsp; k++) {
         XML_Node* s = spDataNodeList[k];
diff --git a/src/thermo/WaterProps.cpp b/src/thermo/WaterProps.cpp
index cc0c6dce1..72ee3540f 100644
--- a/src/thermo/WaterProps.cpp
+++ b/src/thermo/WaterProps.cpp
@@ -287,7 +287,7 @@ doublereal WaterProps::coeffThermalExp_IAPWS(doublereal temp, doublereal press)
     doublereal dens = m_waterIAPWS->density(temp, press, WATER_LIQUID);
     if (dens < 0.0) {
         throw CanteraError("WaterProps::coeffThermalExp_IAPWS",
-                           "Unable to solve for density at T = " + fp2str(temp) + " and P = " + fp2str(press));
+            "Unable to solve for density at T = {} and P = {}", temp, press);
     }
     return m_waterIAPWS->coeffThermExp();
 }
@@ -297,7 +297,7 @@ doublereal WaterProps::isothermalCompressibility_IAPWS(doublereal temp, doublere
     doublereal dens = m_waterIAPWS->density(temp, press, WATER_LIQUID);
     if (dens < 0.0) {
         throw CanteraError("WaterProps::isothermalCompressibility_IAPWS",
-                           "Unable to solve for density at T = " + fp2str(temp) + " and P = " + fp2str(press));
+            "Unable to solve for density at T = {} and P = {}", temp, press);
     }
     return m_waterIAPWS->isothermalCompressibility();
 }
diff --git a/src/thermo/WaterPropsIAPWS.cpp b/src/thermo/WaterPropsIAPWS.cpp
index b6b3a23ac..16aaaaca7 100644
--- a/src/thermo/WaterPropsIAPWS.cpp
+++ b/src/thermo/WaterPropsIAPWS.cpp
@@ -130,7 +130,7 @@ doublereal WaterPropsIAPWS::density(doublereal temperature, doublereal pressure,
                                        "Unstable Branch finder is untested");
                 } else {
                     throw CanteraError("WaterPropsIAPWS::density",
-                                       "unknown state: " + int2str(phase));
+                                       "unknown state: {}", phase);
                 }
             }
         } else {
@@ -188,7 +188,7 @@ doublereal WaterPropsIAPWS::density_const(doublereal pressure,
                                        "Unstable Branch finder is untested");
                 } else {
                     throw CanteraError("WaterPropsIAPWS::density",
-                                       "unknown state: " + int2str(phase));
+                                       "unknown state: {}", phase);
                 }
             }
         } else {
@@ -309,8 +309,8 @@ void WaterPropsIAPWS::corr(doublereal temperature, doublereal pressure,
     densLiq = density(temperature, pressure, WATER_LIQUID, densLiq);
     if (densLiq <= 0.0) {
         throw CanteraError("WaterPropsIAPWS::corr",
-                           "Error occurred trying to find liquid density at (T,P) = "
-                           + fp2str(temperature) + "  " + fp2str(pressure));
+            "Error occurred trying to find liquid density at (T,P) = {}  {}",
+            temperature, pressure);
     }
     setState_TR(temperature, densLiq);
     doublereal gibbsLiqRT = m_phi->gibbs_RT();
@@ -318,8 +318,8 @@ void WaterPropsIAPWS::corr(doublereal temperature, doublereal pressure,
     densGas = density(temperature, pressure, WATER_GAS, densGas);
     if (densGas <= 0.0) {
         throw CanteraError("WaterPropsIAPWS::corr",
-                           "Error occurred trying to find gas density at (T,P) = "
-                           + fp2str(temperature) + "  " + fp2str(pressure));
+            "Error occurred trying to find gas density at (T,P) = {}  {}",
+            temperature, pressure);
     }
     setState_TR(temperature, densGas);
     doublereal gibbsGasRT = m_phi->gibbs_RT();
@@ -333,8 +333,8 @@ void WaterPropsIAPWS::corr1(doublereal temperature, doublereal pressure,
     densLiq = density(temperature, pressure, WATER_LIQUID, densLiq);
     if (densLiq <= 0.0) {
         throw CanteraError("WaterPropsIAPWS::corr1",
-                           "Error occurred trying to find liquid density at (T,P) = "
-                           + fp2str(temperature) + "  " + fp2str(pressure));
+            "Error occurred trying to find liquid density at (T,P) = {}  {}",
+            temperature, pressure);
     }
     setState_TR(temperature, densLiq);
     doublereal prL = m_phi->phiR();
@@ -342,8 +342,8 @@ void WaterPropsIAPWS::corr1(doublereal temperature, doublereal pressure,
     densGas = density(temperature, pressure, WATER_GAS, densGas);
     if (densGas <= 0.0) {
         throw CanteraError("WaterPropsIAPWS::corr1",
-                           "Error occurred trying to find gas density at (T,P) = "
-                           + fp2str(temperature) + "  " + fp2str(pressure));
+            "Error occurred trying to find gas density at (T,P) = {}  {}",
+            temperature, pressure);
     }
     setState_TR(temperature, densGas);
     doublereal prG = m_phi->phiR();
@@ -389,7 +389,7 @@ doublereal WaterPropsIAPWS::psat(doublereal temperature, int waterState)
         setState_TR(temperature, densGas);
     } else {
         throw CanteraError("WaterPropsIAPWS::psat",
-                           "unknown water state input: " + int2str(waterState));
+                           "unknown water state input: {}", waterState);
     }
     return p;
 }
diff --git a/src/thermo/WaterSSTP.cpp b/src/thermo/WaterSSTP.cpp
index a513fa26c..1e321f31a 100644
--- a/src/thermo/WaterSSTP.cpp
+++ b/src/thermo/WaterSSTP.cpp
@@ -383,7 +383,7 @@ doublereal WaterSSTP::dthermalExpansionCoeffdT() const
     doublereal dd = m_sub->density(tt, pres, WATER_LIQUID, dens_save);
     if (dd < 0.0) {
         throw CanteraError("WaterSSTP::dthermalExpansionCoeffdT",
-                           "Unable to solve for the density at T = " + fp2str(tt) + ", P = " + fp2str(pres));
+            "Unable to solve for the density at T = {}, P = {}", tt, pres);
     }
     doublereal vald = m_sub->coeffThermExp();
     m_sub->setState_TR(T, dens_save);
diff --git a/src/transport/HighPressureGasTransport.cpp b/src/transport/HighPressureGasTransport.cpp
index e675f16bd..998e5165e 100644
--- a/src/transport/HighPressureGasTransport.cpp
+++ b/src/transport/HighPressureGasTransport.cpp
@@ -258,7 +258,7 @@ void HighPressureGasTransport::getMultiDiffCoeffs(const size_t ld, doublereal* c
     int ierr = invert(m_Lmatrix, m_nsp);
     if (ierr != 0) {
         throw CanteraError("HighPressureGasTransport::getMultiDiffCoeffs",
-                            string(" invert returned ierr = ")+int2str(ierr));
+                           "invert returned ierr = {}", ierr);
     }
     m_l0000_ok = false; // matrix is overwritten by inverse
     m_lmatrix_soln_ok = false;
diff --git a/src/transport/LiquidTransport.cpp b/src/transport/LiquidTransport.cpp
index 006a7e8b3..94b252191 100644
--- a/src/transport/LiquidTransport.cpp
+++ b/src/transport/LiquidTransport.cpp
@@ -719,7 +719,7 @@ bool LiquidTransport::update_T()
     // Next do a reality check on temperature value
     if (t < 0.0) {
         throw CanteraError("LiquidTransport::update_T()",
-                           "negative temperature "+fp2str(t));
+                           "negative temperature {}", t);
     }
 
     // Compute various direct functions of temperature
diff --git a/src/transport/MixTransport.cpp b/src/transport/MixTransport.cpp
index 211aa814d..69d0f29d5 100644
--- a/src/transport/MixTransport.cpp
+++ b/src/transport/MixTransport.cpp
@@ -128,7 +128,7 @@ void MixTransport::update_T()
     }
     if (t < 0.0) {
         throw CanteraError("MixTransport::update_T",
-                           "negative temperature "+fp2str(t));
+                           "negative temperature {}", t);
     }
     GasTransport::update_T();
     // temperature has changed, so polynomial fits will need to be redone.
diff --git a/src/transport/MultiTransport.cpp b/src/transport/MultiTransport.cpp
index 4f6d79619..162323382 100644
--- a/src/transport/MultiTransport.cpp
+++ b/src/transport/MultiTransport.cpp
@@ -347,12 +347,12 @@ void MultiTransport::getMassFluxes(const doublereal* state1, const doublereal* s
     int info = m_aa.factor();
     if (info) {
         throw CanteraError("MultiTransport::getMassFluxes",
-                           "Error in factorization.  Info = "+int2str(info));
+                           "Error in factorization.  Info = {}", info);
     }
     info = m_aa.solve(fluxes);
     if (info) {
         throw CanteraError("MultiTransport::getMassFluxes",
-                           "Error in linear solve. Info = "+int2str(info));
+                           "Error in linear solve. Info = {}", info);
     }
 
     doublereal pp = pressure_ig();
@@ -403,7 +403,7 @@ void MultiTransport::getMultiDiffCoeffs(const size_t ld, doublereal* const d)
     int ierr = invert(m_Lmatrix, m_nsp);
     if (ierr != 0) {
         throw CanteraError("MultiTransport::getMultiDiffCoeffs",
-                           string(" invert returned ierr = ")+int2str(ierr));
+                           "invert returned ierr = {}", ierr);
     }
     m_l0000_ok = false; // matrix is overwritten by inverse
     m_lmatrix_soln_ok = false;
diff --git a/src/transport/SimpleTransport.cpp b/src/transport/SimpleTransport.cpp
index ffc5a43e5..fe7588589 100644
--- a/src/transport/SimpleTransport.cpp
+++ b/src/transport/SimpleTransport.cpp
@@ -616,7 +616,7 @@ bool SimpleTransport::update_T()
     }
     if (t < 0.0) {
         throw CanteraError("SimpleTransport::update_T",
-                           "negative temperature "+fp2str(t));
+                           "negative temperature {}", t);
     }
 
     // Compute various functions of temperature
diff --git a/src/transport/TransportData.cpp b/src/transport/TransportData.cpp
index d436c7b3a..79d5d7410 100644
--- a/src/transport/TransportData.cpp
+++ b/src/transport/TransportData.cpp
@@ -56,51 +56,49 @@ void GasTransportData::validate(const Species& sp)
     if (geometry == "atom") {
         if (nAtoms != 1) {
             throw CanteraError("GasTransportData::validate",
-                "invalid geometry for species '" + sp.name + "'. 'atom' "
-                "specified, but species contains multiple atoms.");
+                "invalid geometry for species '{}'. 'atom' specified, but "
+                "species contains multiple atoms.", sp.name);
         }
     } else if (geometry == "linear") {
         if (nAtoms == 1) {
             throw CanteraError("GasTransportData::validate",
-                "invalid geometry for species '" + sp.name + "'. 'linear'"
-                " specified, but species only contains one atom.");
+                "invalid geometry for species '{}'. 'linear' specified, but "
+                "species only contains one atom.", sp.name);
         }
     } else if (geometry == "nonlinear") {
         if (nAtoms < 3) {
             throw CanteraError("GasTransportData::validate",
-                "invalid geometry for species '" + sp.name + "'. 'nonlinear'"
-                " specified, but species only contains " + fp2str(nAtoms) +
-                " atoms.");
+                "invalid geometry for species '{}'. 'nonlinear' specified, but "
+                "species only contains {} atoms.", sp.name, nAtoms);
         }
     } else {
         throw CanteraError("GasTransportData::validate",
-                           "invalid geometry for species '" + sp.name + "': '" +
-                           geometry + "'.");
+            "invalid geometry for species '{}': '{}'.", sp.name, geometry);
     }
 
     if (well_depth < 0.0) {
         throw CanteraError("GasTransportData::validate",
-                           "negative well depth for species '" + sp.name + "'.");
+                           "negative well depth for species '{}'.", sp.name);
     }
 
     if (diameter <= 0.0) {
         throw CanteraError("GasTransportData::validate",
-            "negative or zero diameter for species '" + sp.name + "'.");
+            "negative or zero diameter for species '{}'.", sp.name);
     }
 
     if (dipole < 0.0) {
         throw CanteraError("GasTransportData::validate",
-            "negative dipole moment for species '" + sp.name + "'.");
+            "negative dipole moment for species '{}'.", sp.name);
     }
 
     if (polarizability < 0.0) {
         throw CanteraError("GasTransportData::validate",
-            "negative polarizability for species '" + sp.name + "'.");
+            "negative polarizability for species '{}'.", sp.name);
     }
 
     if (rotational_relaxation < 0.0) {
         throw CanteraError("GasTransportData::validate",
-            "negative rotation relaxation number for species '" + sp.name + "'");
+            "negative rotation relaxation number for species '{}'.", sp.name);
     }
 }
 
diff --git a/src/zeroD/Reactor.cpp b/src/zeroD/Reactor.cpp
index b47433ff9..0247478ab 100644
--- a/src/zeroD/Reactor.cpp
+++ b/src/zeroD/Reactor.cpp
@@ -163,12 +163,9 @@ void Reactor::updateState(doublereal* y)
             T -= dT;
             i++;
             if (i > 100) {
-                std::string message = "no convergence";
-                message += "\nU/m = " + fp2str(U / m_mass);
-                message += "\nT = " + fp2str(T);
-                message += "\nrho = " + fp2str(m_mass / m_vol);
-                message += "\n";
-                throw CanteraError("Reactor::updateState", message);
+                throw CanteraError("Reactor::updateState",
+                    "no convergence\nU/m = {}\nT = {}\nrho = {}\n",
+                    U / m_mass, T, m_mass / m_vol);
             }
         }
     } else {
@@ -316,7 +313,7 @@ void Reactor::addSensitivityReaction(size_t rxn)
 {
     if (rxn >= m_kin->nReactions()) {
         throw CanteraError("Reactor::addSensitivityReaction",
-                           "Reaction number out of range ("+int2str(rxn)+")");
+                           "Reaction number out of range ({})", rxn);
     }
 
     network().registerSensitivityReaction(this, rxn,
diff --git a/src/zeroD/Wall.cpp b/src/zeroD/Wall.cpp
index 99f8cd4cd..e11656428 100644
--- a/src/zeroD/Wall.cpp
+++ b/src/zeroD/Wall.cpp
@@ -143,7 +143,7 @@ void Wall::addSensitivityReaction(int leftright, size_t rxn)
 {
     if (rxn >= m_chem[leftright]->nReactions()) {
         throw CanteraError("Wall::addSensitivityReaction",
-                           "Reaction number out of range ("+int2str(rxn)+")");
+                           "Reaction number out of range ({})", rxn);
     }
     if (leftright == 0) {
         m_left->network().registerSensitivityReaction(this, rxn,