fix error in #320

src/methods/property_solvers/Pproperty.jl

@@ -196,6 +196,7 @@ function Pproperty_pure(model,T,prop,z,property::F,rootsolver,phase,abstol,relto
   Psat,vlsat,vvpat = saturation_pressure(model,T,crit = crit)
 
   if !is_unknown(phase)
+
     return __Pproperty(model,T,prop,z,property,rootsolver,phase,abstol,reltol,threaded,Psat)
   end
 
@@ -235,13 +236,13 @@ function __Pproperty(model,T,prop,z,property::F,rootsolver,phase,abstol,reltol,t
     end
     return __Pproperty(model,T,prop,z,property,rootsolver,new_phase,abstol,reltol,threaded,p0)
   end
-  f(p,prop) = property(model,p,T,z,phase = phase,threaded = threaded) - prop
-  prob = Roots.ZeroProblem(f,p0)
+  f(lnp,prop) = property(model,exp(lnp),T,z,phase = phase,threaded = threaded) - prop
+  prob = Roots.ZeroProblem(f,log(p0))
   sol = Roots.solve(prob,rootsolver,p = prop,atol = abstol,rtol = reltol)
-  if !isfinite(sol) || sol < 0
+  if isnan(sol)
     return sol,:failure
   end
-  return sol,phase
+  return exp(sol),phase
 end
 #=
 model = PCSAFT(["propane","dodecane"])

src/methods/property_solvers/multicomponent/flash/general_flash.jl

@@ -386,7 +386,7 @@ function xy_flash_neq(output,model,zbulk,np,input,state::F,μconfig) where F
         i_spec = spec1.k
         output[end-1] = β[i_spec] - val1
     else #general caloric term, valid for enthalpy , entropy, internal energy, gibbs, helmholtz
-        output[end-1] = (val1 - ∑a - p_end*∑v - T*∑s)/val1
+        output[end-1] = (val1 - ∑a - p_end*∑v - T*∑s)/RTinv
     end
 
     if spec2 == temperature
@@ -399,7 +399,7 @@ function xy_flash_neq(output,model,zbulk,np,input,state::F,μconfig) where F
         i_spec = spec2.k
         output[end-1] = β[i_spec] - val2
     else
-        output[end] = (val2 - ∑a - p_end*∑v - T*∑s)/val1
+        output[end] = (val2 - ∑a - p_end*∑v - T*∑s)/RTinv
     end
     return output
 end
@@ -768,16 +768,16 @@ function px_flash_pure(model,p,x,z,spec::F,T0 = nothing) where F
     Ts,vl,vv = saturation_temperature(model,p)
     ∑z = sum(z)
     x1 = SA[1.0]
-    spec_to_vt(model,vl,T,x1,spec)
-    xl = ∑z*spec_to_vt(model,vl,T,x1,spec)
-    xv = ∑z*spec_to_vt(model,vv,T,x1,spec)
+    spec_to_vt(model,vl,Ts,x1,spec)
+    xl = ∑z*spec_to_vt(model,vl,Ts,x1,spec)
+    xv = ∑z*spec_to_vt(model,vv,Ts,x1,spec)
     βv = (x - xl)/(xv - xl)
     if !isfinite(βv)
         return FlashResultInvalid(1,βv)
     elseif βv < 0 || βv > 1
         phase0 = βv < 0 ? :liquid : :vapour
-        T,_phase = _Tproperty(model,p,h/∑z,SA[1.0],spec,T0 = T0,phase = phase0)
-        return FlashResult(model,p,T,∑z,phase = _phase)
+        T,_phase = _Tproperty(model,p,x/∑z,SA[1.0],spec,T0 = T0,phase = phase0)
+        return FlashResult(model,p,T,SA[∑z*one(p)*one(T)],phase = _phase)
     else
         build_flash_result_pure(model,p,Ts,z,vl,vv,βv)
     end
@@ -811,10 +811,10 @@ function tx_flash_pure(model,T,x,z,spec::F,P0 = nothing) where F
         return FlashResultInvalid(1,βv)
     elseif βv < 0 || βv > 1
         phase0 = βv < 0 ? :liquid : :vapour
-        p,_phase = _Pproperty(model,T,x/∑z,SA[1.0],spec,p0 = P0,phase = phase0)
-        return FlashResult(model,p,T,∑z,phase = _phase)
+        p,_phase = _Pproperty(model,T,x/∑z,SA[1.0],spec,p0 = P0)
+        return FlashResult(model,p,T,SA[∑z*one(p)*one(T)],phase = _phase)
     else
-        build_flash_result_pure(model,p,Ts,z,vl,vv,βv)
+        build_flash_result_pure(model,p,T,z,vl,vv,βv)
     end
 end
 

test/test_methods_api.jl

@@ -414,6 +414,13 @@ end
     result = xy_flash(model,spec,z,result0) #perform the flash
     @test Clapeyron.temperature(result) == 200.15
     @test pressure(result) == 101325.0
+
+    #px_flash_pure/tx_flash_pure, 1 phase (#320)
+    model = cPR(["ethane"],idealmodel = ReidIdeal);
+    p = 101325;h = 100; z = Clapeyron.SA[1]; T = Clapeyron.PH.temperature(model,p,h,z)
+    @test enthalpy(model,p,T,z) ≈ h rtol = 1e-6
+    res5 = Clapeyron.tx_flash_pure(model,T,h,z,enthalpy)
+    @test pressure(res5)  ≈ p rtol = 1e-6
 end
 
 @testset "Saturation Methods" begin