simulkade / pvttool

A Matlab toolbox for the PVT calculation using the cubic equations of state (EOS)

PVTtoolbox

A Matlab toolbox for the PVT calculation using the cubic equations of state (EOS)

Getting started

Open Matlab, and run PVTinitialize. Then go to the Examples folder and run the test cases.

pvttool's Issues

Why are the gas phase compressibility factors equal to those of the liquid phase?

Subject: Why are the gas phase compressibility factors equal to those of the liquid phase?

Dear Professor,

I hope this email finds you well. I am writing to express my gratitude for the wonderful code you have developed for me to study. However, while working with the PREOS(mixture1, thermo1) function, I encountered a problem. Specifically, I noticed that the gas phase compressibility factors are equal to those of the liquid phase, even though the component composition results for each phase are correct. The issue arises when I attempt to calculate the density using the compressibility factor.

I am writing to seek your help and advice on this matter. I would greatly appreciate it if you could provide some insights into why the gas phase compressibility factors are equal to those of the liquid phase, and how I can solve the problem of calculating density using the compressibility factor.

Thank you again for your support, and I look forward to hearing from you soon.

Best regards,
RPeter

Stability Testing

Dear Ali,

I am very grateful to study PVTtool. Here I have a question. When I run code of the 'testcase2':

% test case
% liquid-two-phase phase transition
clc; clear;
% Define the components and load pure physsical properties
[component, comp_flag] = addComponents({'CH4', 'C2H6', 'C4H10', ...
'C6H14', 'C10H22', 'C15H32'});
% Define the thermodynamic models
T0 = 300; % [K]
p0 = 100e5; % [Pa]
thermo1 = addThermo();
thermo1.EOS = @PREOS;
mixture1 = addMixture(component, T0, p0);
% Define flash options
options.accuracy = 1e-7;
options.iteration = 100;
% Negative flash
[vapor_y, liquid_x, vapor_frac] = ...
        vleflashnegative(mixture1, thermo1, options)
%mixture1.temperature=T_range(i);
[liquid_z, vapor_z, fugacity, HR] = PREOS(mixture1, thermo1)    

[s1, sl, sv] = stabilityTest(mixture1, thermo1)

I found that the result is:

vapor_y =

    0.8653    0.1097    0.0198    0.0048    0.0004    0.0000


liquid_x =

    0.4053    0.1472    0.1165    0.1114    0.1099    0.1097


vapor_frac =

   -0.5188


liquid_z =

    0.5924


vapor_z =

    0.5924


fugacity =

    1.9419    0.3795    0.0352    0.0037    0.0001    0.0000


HR =

  -3.0503e+04


s1 =

     1     1


sl =

    0.9910


sv =

    1.0000

The problem is

s1 =

     1     1

that means it converged to trivial solution. If so, the system is single-phase but the simulation gives
me two phase data. I think this code did not implement stability analysis to flash calculation directly.
If it is problem, I would be happy to help modify this model.

Thank you so much.

Step by Step

I got a ternary mixture of CH4, N2 and H2 and want to get VLE Diagramms. Could someone explain me how to use the PVT Tool stepbystep. I don't get it. :/
I used Matlab 2019 and tried to first addComponents but the tool only allows one Component and shows me a 0x0 matrix. I also tried it with addMixture('CH4' 'H2' 'N2', 300, 100000) but it didn't work because of to many input arguments.
I would be really grateful if someone could help.

vleflashnegative.m seems cannot be converged.

Hi professor,

Hope this message finds you well.

I try to run your examples, but I find the function vleflashnegative.m seems not work.
For example, the example 'blackoilmodel'：

% test case
% liquid-two-phase phase transition
clc; clear;
% Define the components and load pure physsical properties
z=[0.1, 0.1, 0.1, 0.2, 0.3, 0.2]
[component, comp_flag] = addComponents({'CH4', 'C2H6', 'C4H10', ...
'C6H14', 'C8H18', 'C10H22'});
% Define the thermodynamic models
T = 100+273.15; % [K]
T0= 25+273.15; % [K] surface temperature
p0= 1e5; % [Pa] surface pressure
thermo1 = addThermo();
thermo1.EOS = @preos;
mixture1 = addMixture(component, T0, p0);
mixture1.mole_fraction = z;
% Define flash options
options.accuracy = 1e-7;
options.iteration = 100;
% flash for surface properties
[vapor_y, liquid_x, vapor_frac_surf] = ...
vleflash(mixture1, thermo1, options)
[liquid_z, vapor_z, fugacity, HR] = PREOS(mixture1, thermo1);
R=8.314; % [J/(mol.K)]
mixture1.mole_fraction=liquid_x;
[liquid_z_surf, vapor_z, fugacity, HR] = PREOS(mixture1, thermo1);
V_oil_surface = liquid_z_surf*(1-vapor_frac_surf)RT0/p0;
mixture1.mole_fraction=vapor_y;
[liquid_z, vapor_z_surf, fugacity, HR] = PREOS(mixture1, thermo1);
V_gas_surface = vapor_z_surfvapor_frac_surfR*T0/p0;

n_point=50;
p_range=linspace(p0, 300e5, n_point);
for i=1:n_point
p=p_range(i);
% flash for reservoir properties
mixture1.temperature=T;
mixture1.pressure=p;
mixture1.mole_fraction = z;
[vapor_y, liquid_x, vapor_frac_res] = ...
vleflashnegative(mixture1, thermo1, options) ; %!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!Here! cannot be converged.
if vapor_frac_res<=0
vapor_y=z;
liquid_x=z;
vapor_frac_res=0;
end
mixture1.mole_fraction=liquid_x;
[liquid_z_res, vapor_z, fugacity, HR] = PREOS(mixture1, thermo1);
V_oil_res = liquid_z_res*(1-vapor_frac_res)RT/p;
mixture1.mole_fraction=vapor_y;
[liquid_z, vapor_z_res, fugacity, HR] = PREOS(mixture1, thermo1);
V_gas_res = vapor_z_resvapor_frac_resR*T/p;

%[s1, sl, sv] = stabilityTest(mixture1, thermo1)
Bo(i) = V_oil_res/V_oil_surface
Bg(i) = V_gas_res/V_oil_surface
Rs(i) = (vapor_frac_surf-vapor_frac_res)vapor_z_surfRT0/p0/V_oil_surface
rs1 = (vapor_frac_res-vapor_frac_surf)liquid_z_surfRT0/p0/V_gas_surface;
if rs1<0
rs1=0;
end
rs(i)=rs1;
end

Could you tell me the problem?
AND could you please teach me what is the different between vleflash and vleflashnegative, and their applicable scene?

Asking for calculating bubble point line and dew point line

I am writing to express my sincere gratitude for your excellent work.
As a user, I have found your work to be incredibly helpful in learning phase calculation.
I also have a query regarding the bubble point line and dew point line that I found the code file 'bubblePressure' /
'bubbleTemperature' /'dewPressure'/'dewTemperature' in the path of '\PVTtool-master\Flash' are empty.
Could you please provide the above complete code file and the demo example to calculate and plot bubble point line, dew point line, and isometric line?
Thank you sooooo much!!!