diff --git a/src/NumericalInitialization.cc b/src/NumericalInitialization.cc
index 8e899668297e39120b910e3f6809f80d06078f89..e004adc402caa63b31525903f9ad7eace9fdb4df 100644
--- a/src/NumericalInitialization.cc
+++ b/src/NumericalInitialization.cc
@@ -619,5 +619,8 @@ void
 LoadParamsAndSteadyStateStatement::fillEvalContext(eval_context_t &eval_context) const
 {
   for (const auto & it : content)
-    eval_context[it.first] = stod(it.second);
+    /* We use strtod() instead of stod() because we want overflows and
+       underflows to respectively yield 0 and ±Inf. See also the comment in
+       NumericalConstants.cc */
+    eval_context[it.first] = strtod(it.second.c_str(), nullptr);
 }