diff --git a/include/cantera/base/AnyMap.h b/include/cantera/base/AnyMap.h
index 9b35ad0b5..6bf1da64d 100644
--- a/include/cantera/base/AnyMap.h
+++ b/include/cantera/base/AnyMap.h
@@ -45,6 +45,7 @@ public:
     AnyValue& operator=(AnyValue&& other);
 
     AnyValue& operator[](const std::string& key);
+    const AnyValue& operator[](const std::string& key) const;
 
     bool hasKey(const std::string& key) const;
 
@@ -193,7 +194,8 @@ std::vector<vector_fp>& AnyValue::asVector<vector_fp>(size_t nMin, size_t nMax);
  *     breakfast["waffle"].asDouble();
  * } except (std::exception& err) {
  *     // Exception will be thrown.
- *     // 'breakfast' will have an empty key named "waffle"
+ *     // 'breakfast' will have an empty key named "waffle" unless `breakfast`
+ *     // is a `const AnyMap`.
  * }
  *
  * try {
@@ -232,6 +234,7 @@ public:
     static AnyMap fromYamlString(const std::string& yaml);
 
     AnyValue& operator[](const std::string& key);
+    const AnyValue& operator[](const std::string& key) const;
 
     const AnyValue& at(const std::string& key) const;
 
diff --git a/src/base/AnyMap.cpp b/src/base/AnyMap.cpp
index e0fa4355d..3c6a9b538 100644
--- a/src/base/AnyMap.cpp
+++ b/src/base/AnyMap.cpp
@@ -276,6 +276,11 @@ AnyValue& AnyValue::operator[](const std::string& key)
     return as<AnyMap>()[key];
 }
 
+const AnyValue& AnyValue::operator[](const std::string& key) const
+{
+    return as<AnyMap>().at(key);
+}
+
 bool AnyValue::hasKey(const std::string& key) const {
     return (is<AnyMap>() && as<AnyMap>().hasKey(key));
 }
@@ -367,7 +372,7 @@ std::unordered_map<std::string, const AnyMap*> AnyValue::asMap(
 {
     std::unordered_map<std::string, const AnyMap*> mapped;
     for (const auto& item : asVector<AnyMap>()) {
-        auto key = item.at(name).asString();
+        auto key = item[name].asString();
         if (mapped.count(key)) {
             throw CanteraError("AnyValue::asMap", "Duplicate key '{}'", key);
         }
@@ -400,7 +405,7 @@ void AnyValue::applyUnits(const UnitSystem& units)
             // First item in the list is a units declaration, which applies to
             // the items in the list
             UnitSystem newUnits = units;
-            newUnits.setDefaults(list[0].at("units").asMap<std::string>());
+            newUnits.setDefaults(list[0]["units"].asMap<std::string>());
             list[0].m_data.erase("units");
             for (auto& item : list) {
                 // Any additional units declarations are errors
@@ -560,6 +565,16 @@ AnyValue& AnyMap::operator[](const std::string& key)
     }
 }
 
+const AnyValue& AnyMap::operator[](const std::string& key) const
+{
+    try {
+        return m_data.at(key);
+    } catch (std::out_of_range& err) {
+        throw CanteraError("AnyMap::operator[]",
+            "Key '{}' not found.\nExisting keys: {}", key, keys_str());
+    }
+}
+
 const AnyValue& AnyMap::at(const std::string& key) const
 {
     try {
diff --git a/src/kinetics/Reaction.cpp b/src/kinetics/Reaction.cpp
index 14e95700f..3898d5f4d 100644
--- a/src/kinetics/Reaction.cpp
+++ b/src/kinetics/Reaction.cpp
@@ -373,20 +373,20 @@ void readFalloff(FalloffReaction& R, const XML_Node& rc_node)
 void readFalloff(FalloffReaction& R, const AnyMap& node)
 {
     if (node.hasKey("Troe")) {
-        auto& f = node.at("Troe").as<AnyMap>();
+        auto& f = node["Troe"].as<AnyMap>();
         vector_fp params{
-            f.at("A").asDouble(),
-            f.at("T3").asDouble(),
-            f.at("T1").asDouble(),
+            f["A"].asDouble(),
+            f["T3"].asDouble(),
+            f["T1"].asDouble(),
             f.getDouble("T2", 0.0)
         };
         R.falloff = newFalloff(TROE_FALLOFF, params);
     } else if (node.hasKey("SRI")) {
-        auto& f = node.at("SRI").as<AnyMap>();
+        auto& f = node["SRI"].as<AnyMap>();
         vector_fp params{
-            f.at("A").asDouble(),
-            f.at("B").asDouble(),
-            f.at("C").asDouble(),
+            f["A"].asDouble(),
+            f["B"].asDouble(),
+            f["C"].asDouble(),
             f.getDouble("D", 1.0),
             f.getDouble("E", 0.0)
         };
@@ -411,7 +411,7 @@ void readEfficiencies(ThirdBody& tbody, const AnyMap& node)
 {
     tbody.default_efficiency = node.getDouble("default-efficiency", 1.0);
     if (node.hasKey("efficiencies")) {
-        tbody.efficiencies = node.at("efficiencies").asMap<double>();
+        tbody.efficiencies = node["efficiencies"].asMap<double>();
     }
 }
 
@@ -439,7 +439,7 @@ void setupReaction(Reaction& R, const AnyMap& node)
     // Parse the reaction equation to determine participating species and
     // stoichiometric coefficients
     std::vector<std::string> tokens;
-    tokenizeString(node.at("equation").asString(), tokens);
+    tokenizeString(node["equation"].asString(), tokens);
     tokens.push_back("+"); // makes parsing last species not a special case
 
     size_t last_used = npos; // index of last-used token
@@ -465,7 +465,7 @@ void setupReaction(Reaction& R, const AnyMap& node)
             } else {
                 throw CanteraError("setupReaction", "Error parsing reaction "
                     "string '{}'.\nCurrent token: '{}'\nlast_used: '{}'",
-                    node.at("equation").asString(),
+                    node["equation"].asString(),
                     tokens[i], (last_used == npos) ? "n/a" : tokens[last_used]
                     );
             }
@@ -492,7 +492,7 @@ void setupReaction(Reaction& R, const AnyMap& node)
     // Non-stoichiometric reaction orders
     std::map<std::string, double> orders;
     if (node.hasKey("orders")) {
-        for (const auto& order : node.at("orders").asMap<double>()) {
+        for (const auto& order : node["orders"].asMap<double>()) {
             R.orders[order.first] = order.second;
         }
     }
@@ -531,7 +531,7 @@ void setupElementaryReaction(ElementaryReaction& R, const AnyMap& node,
 {
     setupReaction(R, node);
     R.allow_negative_pre_exponential_factor = node.getBool("negative-A", false);
-    R.rate = readArrhenius(R, node.at("rate-constant"), kin, node.units());
+    R.rate = readArrhenius(R, node["rate-constant"], kin, node.units());
 }
 
 void setupThreeBodyReaction(ThreeBodyReaction& R, const XML_Node& rxn_node)
@@ -547,7 +547,7 @@ void setupThreeBodyReaction(ThreeBodyReaction& R, const AnyMap& node,
     if (R.reactants.count("M") != 1 || R.products.count("M") != 1) {
         throw CanteraError("setupThreeBodyReaction",
             "Reaction equation '{}' does not contain third body 'M'",
-            node.at("equation").asString());
+            node["equation"].asString());
     }
     R.reactants.erase("M");
     R.products.erase("M");
@@ -600,11 +600,11 @@ void setupFalloffReaction(FalloffReaction& R, const AnyMap& node,
     if (third_body == "") {
         throw CanteraError("setupFalloffReaction", "Reactants for reaction "
             "'{}' do not contain a pressure-dependent third body",
-            node.at("equation").asString());
+            node["equation"].asString());
     } else if (R.products.count(third_body) == 0) {
         throw CanteraError("setupFalloffReaction", "Unable to match third body "
             "'{}' in reactants and products of reaction '{}'",
-            third_body, node.at("equation").asString());
+            third_body, node["equation"].asString());
     }
 
     // Remove the dummy species
@@ -619,15 +619,15 @@ void setupFalloffReaction(FalloffReaction& R, const AnyMap& node,
         R.third_body.efficiencies[third_body.substr(2, third_body.size() - 3)] = 1.0;
     }
 
-    if (node.at("type").asString() == "falloff") {
-        R.low_rate = readArrhenius(R, node.at("low-P-rate-constant"), kin,
+    if (node["type"].asString() == "falloff") {
+        R.low_rate = readArrhenius(R, node["low-P-rate-constant"], kin,
                                    node.units(), 1);
-        R.high_rate = readArrhenius(R, node.at("high-P-rate-constant"), kin,
+        R.high_rate = readArrhenius(R, node["high-P-rate-constant"], kin,
                                     node.units());
     } else { // type == "chemically-activated"
-        R.low_rate = readArrhenius(R, node.at("low-P-rate-constant"), kin,
+        R.low_rate = readArrhenius(R, node["low-P-rate-constant"], kin,
                                    node.units());
-        R.high_rate = readArrhenius(R, node.at("high-P-rate-constant"), kin,
+        R.high_rate = readArrhenius(R, node["high-P-rate-constant"], kin,
                                     node.units(), -1);
     }
 
@@ -722,9 +722,9 @@ void setupChebyshevReaction(ChebyshevReaction&R, const AnyMap& node,
     setupReaction(R, node);
     R.reactants.erase("(+M)"); // remove optional third body notation
     R.products.erase("(+M)");
-    const auto& T_range = node.at("temperature-range").asVector<AnyValue>(2);
-    const auto& P_range = node.at("pressure-range").asVector<AnyValue>(2);
-    auto& vcoeffs = node.at("data").asVector<vector_fp>();
+    const auto& T_range = node["temperature-range"].asVector<AnyValue>(2);
+    const auto& P_range = node["pressure-range"].asVector<AnyValue>(2);
+    auto& vcoeffs = node["data"].asVector<vector_fp>();
     Array2D coeffs(vcoeffs.size(), vcoeffs[0].size());
     for (size_t i = 0; i < coeffs.nRows(); i++) {
         for (size_t j = 0; j < coeffs.nColumns(); j++) {
@@ -923,7 +923,7 @@ unique_ptr<Reaction> newReaction(const AnyMap& node, const Kinetics& kin)
 {
     std::string type = "elementary";
     if (node.hasKey("type")) {
-        type = node.at("type").asString();
+        type = node["type"].asString();
     }
 
     if (type == "elementary") {
diff --git a/src/thermo/Species.cpp b/src/thermo/Species.cpp
index 81240a095..cb7487ca0 100644
--- a/src/thermo/Species.cpp
+++ b/src/thermo/Species.cpp
@@ -84,11 +84,11 @@ shared_ptr<Species> newSpecies(const XML_Node& species_node)
 
 unique_ptr<Species> newSpecies(const AnyMap& node)
 {
-    unique_ptr<Species> s(new Species(node.at("name").asString(),
-                                      node.at("composition").asMap<double>()));
+    unique_ptr<Species> s(new Species(node["name"].asString(),
+                                      node["composition"].asMap<double>()));
 
     if (node.hasKey("thermo")) {
-        s->thermo = newSpeciesThermo(node.at("thermo").as<AnyMap>());
+        s->thermo = newSpeciesThermo(node["thermo"].as<AnyMap>());
     } else {
         s->thermo.reset(new SpeciesThermoInterpType());
     }
diff --git a/src/thermo/SpeciesThermoFactory.cpp b/src/thermo/SpeciesThermoFactory.cpp
index ea4dc965d..33808a183 100644
--- a/src/thermo/SpeciesThermoFactory.cpp
+++ b/src/thermo/SpeciesThermoFactory.cpp
@@ -153,7 +153,7 @@ void setupNasaPoly(NasaPoly2& thermo, const AnyMap& node)
 {
     setupSpeciesThermo(thermo, node);
     vector_fp Tranges = node.convertVector("temperature-ranges", "K", 2, 3);
-    const auto& data = node.at("data").asVector<vector_fp>(Tranges.size()-1);
+    const auto& data = node["data"].asVector<vector_fp>(Tranges.size()-1);
     for (const auto& poly : data) {
         if (poly.size() != 7) {
             throw CanteraError("setupNasaPoly", "Wrong number of coefficients "
@@ -249,7 +249,7 @@ void setupShomatePoly(ShomatePoly2& thermo, const AnyMap& node)
 {
     setupSpeciesThermo(thermo, node);
     vector_fp Tranges = node.convertVector("temperature-ranges", "K", 2, 3);
-    const auto& data = node.at("data").asVector<vector_fp>(Tranges.size()-1);
+    const auto& data = node["data"].asVector<vector_fp>(Tranges.size()-1);
     for (const auto& poly : data) {
         if (poly.size() != 7) {
             throw CanteraError("setupShomatePoly", "Wrong number of coefficients "
@@ -349,7 +349,7 @@ void setupNasa9Poly(Nasa9PolyMultiTempRegion& thermo, const AnyMap& node)
 {
     setupSpeciesThermo(thermo, node);
     vector_fp Tranges = node.convertVector("temperature-ranges", "K", 2, 999);
-    const auto& data = node.at("data").asVector<vector_fp>(Tranges.size()-1);
+    const auto& data = node["data"].asVector<vector_fp>(Tranges.size()-1);
     map<double, vector_fp> regions;
     for (size_t i = 0; i < data.size(); i++) {
         if (data[i].size() != 9) {
@@ -370,7 +370,7 @@ void setupMu0(Mu0Poly& thermo, const AnyMap& node)
     bool dimensionless = node.getBool("dimensionless", false);
     double h0 = node.convertMolarEnergy("h0", "J/kmol", 0.0);
     map<double, double> T_mu;
-    for (const auto& item : node.at("data")) {
+    for (const auto& item : node["data"]) {
         double T = node.units().convert(fpValueCheck(item.first), "K");
         if (dimensionless) {
             T_mu[T] = item.second.asDouble() * GasConstant * T;
@@ -438,7 +438,7 @@ SpeciesThermoInterpType* newSpeciesThermoInterpType(const XML_Node& thermo)
 
 unique_ptr<SpeciesThermoInterpType> newSpeciesThermo(const AnyMap& node)
 {
-    std::string model = node.at("model").asString();
+    std::string model = node["model"].asString();
     if (model == "NASA7") {
         unique_ptr<NasaPoly2> thermo(new NasaPoly2());
         setupNasaPoly(*thermo, node);