Add VLE method for pures

See #2

include/teqp/algorithms/VLE.hpp

@@ -61,3 +61,26 @@ public:
         return J;
     }
 };
+
+template<typename Residual, typename Scalar>
+auto do_pure_VLE(Residual &resid, Scalar T, Scalar rhoL, Scalar rhoV) {
+    auto rhovec = (Eigen::ArrayXd(2) << rhoL, rhoV).finished();
+    auto r0 = resid.call(rhovec);
+    auto J = resid.Jacobian(rhovec);
+    for (int iter = 0; iter < 100; ++iter){
+        if (iter > 0) {
+            r0 = resid.call(rhovec);
+            J = resid.Jacobian(rhovec); 
+        }
+        auto v = J.matrix().colPivHouseholderQr().solve(-r0.matrix()).array().eval();
+        auto rhovecnew = (rhovec + v).eval();
+        
+        // If the solution has stopped improving, stop. The change in rhovec is equal to v in infinite precision, but 
+        // not when finite precision is involved, so use the minimum non-denormal float as the determination of whether
+        // the values are done changing
+        if (((rhovecnew - rhovec).cwiseAbs() < std::numeric_limits<Scalar>::min()).all()) {
+            break;
+        }
+        rhovec += v;
+    }
+}
\ No newline at end of file

src/tests/catch_tests.cxx

@@ -242,4 +242,6 @@ TEST_CASE("Test pure VLE", "") {
     auto rhovec1 = rhovec + v;
     auto r1 = resid.call(rhovec1);
     CHECK((r0.cwiseAbs() > r1.cwiseAbs()).eval().all());
+
+    do_pure_VLE(resid, T, 22834.056386882046, 1025.106554560764);
 }
\ No newline at end of file