%build a Pxy diagram for to substances at equilibrium using Rault's law and
%Antoine's equation
T=input('enter the temperature constant of the system: ','s');
A1=input('enter the A coefficient for the first compound: ','s');
B1=input('enter the B coefficient for the first compound: ','s');
C1=input('enter the C coefficient for the first compound: ','s');
A2=input('enter the A coefficient for the second compound: ','s');
B2=input('enter the B coefficient for the second compound: ','s');
C2=input('enter the C coefficient for the second compound: ','s');
x1=input('enter the concentration of the second compound when it begins to boil: ','s');
y1=input('enter the concentration of the first compound when it begins to condense: ','s');
%Pxy bubble uses the given values for a binary two phase system as well as
%the constant pressure and the concentration of  initial dew point.
%   T=temperature. 
%A,B&C, 1&2 are all constant coefficients for the involved compounds 1&2.
%x1 & y1 are the initial concentrations of the compounds.
P1sat= exp(A1-B1/(T+C1));
P2sat= exp(A2-B2/(T+C2));
m= input('would you like to calculate the bubble or dew point?: ','s');
        
        if strcmp(m,'bubble');%if bubble is entered the if statement will solve for the bubble point
             
                P=x1*P1sat+P2sat-x1*P2sat;%antoine's equation
                y2=1-(P1sat*x1)/P;%using Raoult's law
                y1=(P1sat*x1)/P;
                fprintf('the pressure at the bubble point is: ,%s\n',P)%display results
                fprintf('the concentration of y1 at the bubble point is: ,%s\n',y2)
                    %display results
              
         else
             P= P1sat/y1+P2sat/(1-y1);%antoine's equation
             x2=1-(y1*P)/P1sat;%using raoult's law
             x1=(y1*P)/P1sat;
             fprintf('the pressure at the bubble point is: ,%s\n',P)%display results
             fprintf('the concentration of y1 at the bubble point is: ,%s\n',x2)
        end
         i=1;
                for x=0:.1:1;%build a plot from P2sat to P1sat
                    y=0:.1:1;
                    x(i)=(P-P2sat)/(P1sat-P2sat);
                    y(i)=(x1*P1sat/P);
                    i=i+1;
                end
           x=0:.1:1;%set length
           y=0:.1:1;
           figure(2)
           plot(x,P,y,P)

