diff --git a/include/cantera/thermo/PureFluidPhase.h b/include/cantera/thermo/PureFluidPhase.h
index fadb13e03..000c439d0 100644
--- a/include/cantera/thermo/PureFluidPhase.h
+++ b/include/cantera/thermo/PureFluidPhase.h
@@ -161,6 +161,10 @@ public:
     virtual std::string report(bool show_thermo=true,
                                doublereal threshold=1e-14) const;
 
+    virtual bool compatibleWithMultiPhase() const {
+        return false;
+    }
+
 protected:
     //! Main call to the tpx level to set the state of the system
     /*!
diff --git a/include/cantera/thermo/ThermoPhase.h b/include/cantera/thermo/ThermoPhase.h
index fdfffbd7c..d77c77806 100644
--- a/include/cantera/thermo/ThermoPhase.h
+++ b/include/cantera/thermo/ThermoPhase.h
@@ -1290,6 +1290,13 @@ public:
      */
     bool getElementPotentials(doublereal* lambda) const;
 
+    //! Indicates whether this phase type can be used with class MultiPhase for
+    //! equilibrium calculations. Returns `false` for special phase types which
+    //! already represent multi-phase mixtures, namely PureFluidPhase.
+    virtual bool compatibleWithMultiPhase() const {
+        return true;
+    }
+
     //@}
     /// @name Critical State Properties.
     /// These methods are only implemented by subclasses that implement
diff --git a/include/cantera/thermo/WaterSSTP.h b/include/cantera/thermo/WaterSSTP.h
index 46e6a32ff..4a239e97d 100644
--- a/include/cantera/thermo/WaterSSTP.h
+++ b/include/cantera/thermo/WaterSSTP.h
@@ -192,6 +192,10 @@ public:
 
     virtual doublereal satPressure(doublereal t);
 
+    virtual bool compatibleWithMultiPhase() const {
+        return false;
+    }
+
     //! Return the fraction of vapor at the current conditions
     /*!
      * Below Tcrit, this routine will always return 0, by definition of the
diff --git a/interfaces/cython/cantera/test/test_mixture.py b/interfaces/cython/cantera/test/test_mixture.py
index 2f0afa11c..dbeab3981 100644
--- a/interfaces/cython/cantera/test/test_mixture.py
+++ b/interfaces/cython/cantera/test/test_mixture.py
@@ -187,3 +187,8 @@ class TestMixture(utilities.CanteraTest):
             x.foobar = 300
         with self.assertRaises(AttributeError):
             x.foobar
+
+    def test_invalid_phase_type(self):
+        water = ct.Water()
+        with self.assertRaises(Exception):
+            self.mix = ct.Mixture([(self.phase1, 1.0), (water, 2.0)])
diff --git a/src/equil/MultiPhase.cpp b/src/equil/MultiPhase.cpp
index 456616f65..64ee08bc6 100644
--- a/src/equil/MultiPhase.cpp
+++ b/src/equil/MultiPhase.cpp
@@ -89,6 +89,11 @@ void MultiPhase::addPhase(ThermoPhase* p, doublereal moles)
                            "phases cannot be added after init() has been called.");
     }
 
+    if (!p->compatibleWithMultiPhase()) {
+        throw CanteraError("MultiPhase::addPhase", "Phase '{}'' is not "
+            "compatible with MultiPhase equilibrium solver", p->name());
+    }
+
     // save the pointer to the phase object
     m_phase.push_back(p);