#include <stdio.h>
#include <math.h>
#include <conio.h>
#define R 8.314
int main()
{
	int comp,i,j,k,a,b;
	double tempc,tempk,pressdp,pressbp,vi,vio,est,f,fo,df,**mtrx,pressorg,**invv,**tcrit,**zcrit,**wcrit,**vcrit,*gammab,tempb,tempd,*gammad,**pcrit,*phid,*phib,**bij,**delta,*x,*y,*ki,*z,*xo,*yo,*phi,pjsat,tempco,pjsato,*tsat,*gamma,*psat,frac=0.0,press=0.0,sum1=0.0,sum2=0.0,sum3=0.0,presso,*gammao,b0,b1;

	FILE *input;
	input=fopen("data.txt","r");
	fscanf(input,"%d\n",&comp);

	x=(double*)malloc(comp*sizeof(double));
	y=(double*)malloc(comp*sizeof(double));
	phi=(double*)malloc(comp*sizeof(double));
	gamma=(double*)malloc(comp*sizeof(double));
	psat=(double*)malloc(comp*sizeof(double));
	tsat=(double*)malloc(comp*sizeof(double));

	mtrx=(double**)malloc(comp * sizeof (double *));
	for(i=0;i<comp;i++)
		mtrx[i]=(double *) malloc ((9+comp) * sizeof (double));
	
		for(i=0;i<comp;i++)
	{
		for(j=0;j<9+comp;j++)
			fscanf(input,"%lf",&mtrx[i][j]);
	}

	invv=(double**)malloc(comp * sizeof (double *));
	for(i=0;i<comp;i++)
		invv[i]=(double *) malloc ((comp) * sizeof (double));

	tcrit=(double**)malloc(comp * sizeof (double *));
	for(i=0;i<comp;i++)
		tcrit[i]=(double *) malloc ((comp) * sizeof (double));

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		tcrit[i][j]=pow((mtrx[i][1]*mtrx[j][1]),0.5);
	}

	zcrit=(double**)malloc(comp * sizeof (double *));
	for(i=0;i<comp;i++)
		zcrit[i]=(double *) malloc ((comp) * sizeof (double));

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		zcrit[i][j]=(mtrx[i][0]+mtrx[j][0])/2.0;
	}

	wcrit=(double**)malloc(comp * sizeof (double *));
	for(i=0;i<comp;i++)
			wcrit[i]=(double *) malloc ((comp) * sizeof (double));
	
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		wcrit[i][j]=(mtrx[i][3]+mtrx[j][3])/2.0;
	}

	vcrit=(double**)malloc(comp * sizeof (double *));
	for(i=0;i<comp;i++)
		vcrit[i]=(double *) malloc ((comp) * sizeof (double));
	
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		vcrit[i][j]=pow(((pow(mtrx[i][4],(1.0/3.0))+pow(mtrx[j][4],(1.0/3.0)))/2.0),3.0);
	}

	pcrit=(double**)malloc(comp * sizeof (double *));
	for(i=0;i<comp;i++)
		pcrit[i]=(double *) malloc ((comp) * sizeof (double));
	
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		pcrit[i][j]=(zcrit[i][j]*R*1000.0*tcrit[i][j])/vcrit[i][j];
	}

	bij=(double**)malloc(comp * sizeof (double *));
	for(i=0;i<comp;i++)
		bij[i]=(double *) malloc ((comp) * sizeof (double));
	
	delta=(double**)malloc(comp * sizeof (double *));
	for(i=0;i<comp;i++)
		delta[i]=(double *) malloc ((comp) * sizeof (double));


	printf("*************************************************************************\n");
	printf("*			Bogazici University				*\n");
	printf("*		Department of Chemical Engineering			*\n");
	printf("*	     ChE 321 - Introduction to Thermodynamics			*\n");
	printf("*	             Fall 2005  -  Term Project				*\n*									*\n");
	printf("*	  prepared by Ali Riza Tanriverdi - Murat Sezgin		*\n");
	printf("*************************************************************************\n\n");
	printf("This program calculates Dew/Bubble Temperatures and Dew/Bubble Pressures for the given T/P and x/y values in addition the user can perform a flash calculation. The program can calculate these values by using three types of equations: Raoult's Law, Modified Raoult's Law and Gamma/Fi Calculation(the real case).\n");


	printf("\n	Which computation do you want to make ?\n");
	printf("	(1) Bubble Pressure Calculation\n	(2) Dew Pressure Calculation\n	(3) Bubble Temperature Calculation\n	(4) Dew Temperature Calculation\n	(5) P,T Flash Calculation\n\n");
	printf("	Please, enter your choice:");
	scanf("%d",&a);

	
	printf("\n	Which correalation formula do you want to use ?\n");
	printf("	(1) Raoult's Law\n	(2) Modified Raoult's Law with Wilson Equation\n	(3) Gamma-Fi Formulation\n\n");
	printf("	Please, enter your choice:");
	scanf("%d",&b);

	switch(a)
	{
	case 1:
		switch(b)
		{
		case 3:


//BUBBLE P
	printf("Please, enter temperature(in Celcius):");
	scanf("%lf",&tempc);
	tempk=tempc+273.15;

	for(i=0;i<comp-1;i++)
	{
		printf("Please, enter x%d:",i+1);
		scanf("%lf",&x[i]);
		frac=x[i]+frac;
	}
	
	x[comp-1]=1.0-frac;
	
	for(i=0;i<comp;i++)
		phi[i]=1.0;

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];

		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
			sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}

	for(i=0;i<comp;i++)
		press+=x[i]*gamma[i]*psat[i]/phi[i];  


while(1)
	{
	for(i=0;i<comp;i++)
		y[i]=x[i]*gamma[i]*psat[i]/(phi[i]*press);

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)	
		{	
			b0=0.083-(0.422/(pow((tempk/tcrit[i][j]),(1.6))));
			b1=0.139-(0.172/(pow((tempk/tcrit[i][j]),(4.2))));
			bij[i][j]=((R*1000.0*tcrit[i][j])/pcrit[i][j])*(b0+wcrit[i][j]*b1);
		}
	}

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			delta[j][i]=2.0*bij[j][i]-bij[j][j]-bij[i][i];
	}
	sum1=0.0;

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{
			for(k=0;k<comp;k++)
				sum1+=y[j]*y[k]*(2.0*delta[j][i]-delta[j][k]);
		}
		phi[i]=exp((bij[i][i]*(press-psat[i])+0.5*press*sum1)/(R*1000.0*tempk));
		sum1=0.0;
	}

	presso=press;
	press=0.0;
	for(i=0;i<comp;i++)
		press+=x[i]*gamma[i]*psat[i]/phi[i];  

	if(press-presso<0.1&&presso-press>-0.1)
	break;
	}
		printf("Pressure:%f kPa\n",press);
		for(i=0;i<comp;i++)
			printf("y%d:%f   gamma%d:%f   fi%d:%f\n",i+1,y[i],i+1,gamma[i],i+1,phi[i]);

		break;
		

 		case 2:
			printf("Please, enter temperature(in Celcius):");
	scanf("%lf",&tempc);
	tempk=tempc+273.15;

	for(i=0;i<comp-1;i++)
	{
		printf("Please, enter x%d:",i+1);
		scanf("%lf",&x[i]);
		frac=x[i]+frac;
	}
	
	x[comp-1]=1.0-frac;
	
	for(i=0;i<comp;i++)
		phi[i]=1.0;

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];

		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
			sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}

	for(i=0;i<comp;i++)
		press+=x[i]*gamma[i]*psat[i]/phi[i];  


while(1)
	{
	for(i=0;i<comp;i++)
		y[i]=x[i]*gamma[i]*psat[i]/(phi[i]*press);

	presso=press;
	press=0.0;
	for(i=0;i<comp;i++)
		press+=x[i]*gamma[i]*psat[i]/phi[i];  

	if(press-presso<0.1&&presso-press>-0.1)
	break;
	}
		printf("Pressure:%f kPa\n",press);
		for(i=0;i<comp;i++)
			printf("y%d:%f   gamma%d:%f   fi%d:%f\n",i+1,y[i],i+1,gamma[i],i+1,phi[i]);

break;


	case 1:
			printf("Please, enter temperature(in Celcius):");
	scanf("%lf",&tempc);
	tempk=tempc+273.15;

	for(i=0;i<comp-1;i++)
	{
		printf("Please, enter x%d:",i+1);
		scanf("%lf",&x[i]);
		frac=x[i]+frac;
	}
	
	x[comp-1]=1.0-frac;
	
	for(i=0;i<comp;i++)
		phi[i]=1.0;
	for(i=0;i<comp;i++)
		gamma[i]=1.0;

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	for(i=0;i<comp;i++)
		press+=x[i]*gamma[i]*psat[i]/phi[i];  


while(1)
	{
	for(i=0;i<comp;i++)
		y[i]=x[i]*gamma[i]*psat[i]/(phi[i]*press);

	presso=press;
	press=0.0;
	for(i=0;i<comp;i++)
		press+=x[i]*gamma[i]*psat[i]/phi[i];  

	if(press-presso<0.1&&presso-press>-0.1)
	break;
	}
		printf("Pressure:%f kPa\n",press);
		for(i=0;i<comp;i++)
			printf("y%d:%f   gamma%d:%f   fi%d:%f\n",i+1,y[i],i+1,gamma[i],i+1,phi[i]);

break;
		}

		break;

case 2:
	switch(b)
	{
	case 3:



//  DEW P 

	printf("Please, enter temperature(in Celcius):");
	scanf("%lf",&tempc);
	tempk=tempc+273.15;

	frac=0.0;
	for(i=0;i<comp-1;i++)
	{
		printf("Please, enter y%d:",i+1);
		scanf("%lf",&y[i]);
		frac=y[i]+frac;
	}
	y[i]=1.0-frac;
	for(i=0;i<comp;i++)
		phi[i]=1.0;
	for(i=0;i<comp;i++)
		gamma[i]=1.0;
	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	press=0.0;
	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));

	press=1.0/press;

	for(i=0;i<comp;i++)
			x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
	
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;		
	}
	press=0.0;

	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));
	
	press=1.0/press;

while(1)
{
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{	
			b0=0.083-(0.422/(pow((tempk/tcrit[i][j]),(1.6))));
			b1=0.139-(0.172/(pow((tempk/tcrit[i][j]),(4.2))));
			bij[i][j]=((R*1000.0*tcrit[i][j])/pcrit[i][j])*(b0+wcrit[i][j]*b1);    
		}
	}

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			delta[j][i]=2.0*bij[j][i]-bij[j][j]-bij[i][i];
	}

	sum1=0.0;

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{
			for(k=0;k<comp;k++)
				sum1+=y[j]*y[k]*(2.0*delta[j][i]-delta[j][k]);
		}
		
		phi[i]=exp((bij[i][i]*(press-psat[i])+0.5*press*sum1)/(R*1000.0*tempk));
		sum1=0.0;
	}

	while(1)
	{
		for(i=0;i<comp;i++)
			x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

		gammao=(double*)malloc(comp*sizeof(double));
	for(i=0;i<comp;i++)
		gammao[i]=gamma[i];

		sum1=0.0;
		for(i=0;i<comp;i++)
			sum1+=x[i];

		for(i=0;i<comp;i++)
			x[i]=x[i]/sum1;
		
		for(i=0;i<comp;i++)
		{
			for(j=0;j<comp;j++)
				invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
		}
		sum1=0.0;

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
	
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
	
			sum2+=(x[k]*invv[k][i])/sum3;
		}

		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}
	
	j=0;
	for(i=0;i<comp;i++)
	{
		if(gammao[i]-gamma[i]<0.1&&gammao[i]-gamma[i]>-0.1)
		j++;
	}
	if(j==comp)
		break;
	}
	
	presso=press;
	press=0.0;
	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));
	
	press=1.0/press;
	
	if(press-presso<1.0&&presso-press>-1.0)
	break;
	}
	
	printf("Pressure is %f kPa\n",press);
	for(i=0;i<comp;i++)
		printf("x%d:%f   gamma%d:%f   fi%d:%f\n",i+1,x[i],i+1,gamma[i],i+1,phi[i]);
	break;
	

	case 2:

		printf("Please, enter temperature(in Celcius):");
	scanf("%lf",&tempc);
	tempk=tempc+273.15;

	frac=0.0;
	for(i=0;i<comp-1;i++)
	{
		printf("Please, enter y%d:",i+1);
		scanf("%lf",&y[i]);
		frac=y[i]+frac;
	}
	y[i]=1.0-frac;
	for(i=0;i<comp;i++)
		phi[i]=1.0;
	for(i=0;i<comp;i++)
		gamma[i]=1.0;
	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	press=0.0;
	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));

	press=1.0/press;

	for(i=0;i<comp;i++)
			x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
	
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;		
	}
	press=0.0;

	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));
	
	press=1.0/press;

while(1)
{
	while(1)
	{
		for(i=0;i<comp;i++)
			x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

		gammao=(double*)malloc(comp*sizeof(double));
	for(i=0;i<comp;i++)
		gammao[i]=gamma[i];

		sum1=0.0;
		for(i=0;i<comp;i++)
			sum1+=x[i];

		for(i=0;i<comp;i++)
			x[i]=x[i]/sum1;
		
		for(i=0;i<comp;i++)
		{
			for(j=0;j<comp;j++)
				invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
		}
		sum1=0.0;

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
	
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
	
			sum2+=(x[k]*invv[k][i])/sum3;
		}

		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}
	
	j=0;
	for(i=0;i<comp;i++)
	{
		if(gammao[i]-gamma[i]<0.1&&gammao[i]-gamma[i]>-0.1)
		j++;
	}
	if(j==comp)
		break;
	}
	
	presso=press;
	press=0.0;
	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));
	
	press=1.0/press;
	
	if(press-presso<1.0&&presso-press>-1.0)
	break;
	}
	
	printf("Pressure is %f kPa\n",press);
	for(i=0;i<comp;i++)
		printf("x%d:%f   gamma%d:%f   fi%d:%f\n",i+1,x[i],i+1,gamma[i],i+1,phi[i]);
	break;




	case 1:
		printf("Please, enter temperature(in Celcius):");
	scanf("%lf",&tempc);
	tempk=tempc+273.15;

	frac=0.0;
	for(i=0;i<comp-1;i++)
	{
		printf("Please, enter y%d:",i+1);
		scanf("%lf",&y[i]);
		frac=y[i]+frac;
	}
	y[i]=1.0-frac;
	for(i=0;i<comp;i++)
		phi[i]=1.0;
	for(i=0;i<comp;i++)
		gamma[i]=1.0;
	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	press=0.0;
	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));

	press=1.0/press;

	for(i=0;i<comp;i++)
			x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

	press=0.0;

	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));
	
	press=1.0/press;

while(1)
{
	while(1)
	{
		for(i=0;i<comp;i++)
			x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

		gammao=(double*)malloc(comp*sizeof(double));
	for(i=0;i<comp;i++)
		gammao[i]=gamma[i];

		sum1=0.0;
		for(i=0;i<comp;i++)
			sum1+=x[i];

		for(i=0;i<comp;i++)
			x[i]=x[i]/sum1;
			
	j=0;
	for(i=0;i<comp;i++)
	{
		if(gammao[i]-gamma[i]<0.1&&gammao[i]-gamma[i]>-0.1)
		j++;
	}
	if(j==comp)
		break;
	}
	
	presso=press;
	press=0.0;
	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));
	
	press=1.0/press;
	
	if(press-presso<1.0&&presso-press>-1.0)
	break;
	}
	
	printf("Pressure is %f kPa\n",press);
	for(i=0;i<comp;i++)
		printf("x%d:%f   gamma%d:%f   fi%d:%f\n",i+1,x[i],i+1,gamma[i],i+1,phi[i]);
	break;
	}
	break;

	case 3:
		switch(b)
		{
		case 3:
	
//BUBBLE T 

	printf("Please, enter pressure(in kPa):");
	scanf("%lf",&press);

	for(i=0;i<comp-1;i++)
	{
		printf("Please, enter x%d:",i+1);
		scanf("%lf",&x[i]);
		frac=x[i]+frac;
	}
	
	x[i]=1.0-frac;
	for(i=0;i<comp;i++)
		phi[i]=1.0;

	for(i=0;i<comp;i++)
		tsat[i]=(-1.0*mtrx[i][8]+(mtrx[i][7]/(-1.0*log(press)+mtrx[i][6])));
	
	tempc=0.0;
	for(i=0;i<comp;i++)
		tempc+=x[i]*tsat[i];
	
	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));
	
	tempk=tempc+273.15;


	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}
	
	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}


	
	pjsat=psat[0];//the first entry of the data file is selected for the calculation of Pjsat
	sum1=0.0;
	pjsato=psat[0];

	for(i=0;i<comp;i++)
		sum1+=(x[i]*gamma[i]/phi[i])*(psat[i]/pjsat);
	
	pjsat=press/sum1;
		
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

while(1)
{
	
	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	for(i=0;i<comp;i++)
		y[i]=x[i]*gamma[i]*psat[i]/(phi[i]*press);

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{	
			b0=0.083-(0.422/(pow((tempk/tcrit[i][j]),(1.6))));
			b1=0.139-(0.172/(pow((tempk/tcrit[i][j]),(4.2))));
			bij[i][j]=((R*1000.0*tcrit[i][j])/pcrit[i][j])*(b0+wcrit[i][j]*b1);    
		}
	}


	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			delta[j][i]=2.0*bij[j][i]-bij[j][j]-bij[i][i];
	}
	sum1=0.0;

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{
			for(k=0;k<comp;k++)
				sum1+=y[j]*y[k]*(2.0*delta[j][i]-delta[j][k]);
		}
		phi[i]=exp((bij[i][i]*(press-psat[i])+0.5*press*sum1)/(R*1000.0*tempk));
		sum1=0.0;
	}

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
	
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}
	for(i=0;i<comp;i++)
	{
		sum1+=(x[i]*gamma[i]/phi[i])*(psat[i]/pjsat);
	}
		pjsat=press/sum1;


	tempco=tempc;
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

	if(tempco-tempc<0.1&&tempco-tempc>-0.1)
	break;
	}

	printf("Temperature is %f\n",tempc);
	for(i=0;i<comp;i++)
		printf("y%d:%f   gamma%d:%f   fi%d:%f\n",i+1,y[i],i+1,gamma[i],i+1,phi[i]);
	
break;

	case 2:
			printf("Please, enter pressure(in kPa):");
	scanf("%lf",&press);

	for(i=0;i<comp-1;i++)
	{
		printf("Please, enter x%d:",i+1);
		scanf("%lf",&x[i]);
		frac=x[i]+frac;
	}
	
	x[i]=1.0-frac;
	for(i=0;i<comp;i++)
		phi[i]=1.0;

	for(i=0;i<comp;i++)
		tsat[i]=(-1.0*mtrx[i][8]+(mtrx[i][7]/(-1.0*log(press)+mtrx[i][6])));
	
	tempc=0.0;
	for(i=0;i<comp;i++)
		tempc+=x[i]*tsat[i];
	
	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));
	
	tempk=tempc+273.15;


	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}
	
	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}


	
	pjsat=psat[0];//the first entry of the data file is selected for the calculation of Pjsat
	sum1=0.0;
	pjsato=psat[0];

	for(i=0;i<comp;i++)
		sum1+=(x[i]*gamma[i]/phi[i])*(psat[i]/pjsat);
	
	pjsat=press/sum1;
		
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

while(1)
{
	
	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	for(i=0;i<comp;i++)
		y[i]=x[i]*gamma[i]*psat[i]/(phi[i]*press);


	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
	
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}

	for(i=0;i<comp;i++)
	{
		sum1+=(x[i]*gamma[i]/phi[i])*(psat[i]/pjsat);
	}
		pjsat=press/sum1;


	tempco=tempc;
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

	if(tempco-tempc<0.1&&tempco-tempc>-0.1)
	break;
	}

	printf("Temperature is %f\n",tempc);
	for(i=0;i<comp;i++)
		printf("y%d:%f   gamma%d:%f   fi%d:%f\n",i+1,y[i],i+1,gamma[i],i+1,phi[i]);


break;
	case 1:
			printf("Please, enter pressure(in kPa):");
	scanf("%lf",&press);

	for(i=0;i<comp-1;i++)
	{
		printf("Please, enter x%d:",i+1);
		scanf("%lf",&x[i]);
		frac=x[i]+frac;
	}
	
	x[i]=1.0-frac;
	for(i=0;i<comp;i++)
	{
		phi[i]=1.0;
		gamma[i]=1.0;
	}

	
	for(i=0;i<comp;i++)
		tsat[i]=(-1.0*mtrx[i][8]+(mtrx[i][7]/(-1.0*log(press)+mtrx[i][6])));	
	tempc=0.0;
	for(i=0;i<comp;i++)
		tempc+=x[i]*tsat[i];

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	tempk=tempc+273.15;


	pjsat=psat[0];//the first entry of the data file is selected for the calculation of Pjsat
	sum1=0.0;
	pjsato=psat[0];

	for(i=0;i<comp;i++)
		sum1+=(x[i]*gamma[i]/phi[i])*(psat[i]/pjsat);
	
	pjsat=press/sum1;
		
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;


while(1)
{
	
	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));


	for(i=0;i<comp;i++)
		y[i]=x[i]*gamma[i]*psat[i]/(phi[i]*press);




	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		sum1+=(x[i]*gamma[i]/phi[i])*(psat[i]/pjsat);
	}
		pjsat=press/sum1;


	tempco=tempc;
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;


	if(tempco-tempc<0.1&&tempco-tempc>-0.1)
	break;
	}

	printf("Temperature is %f\n",tempc);
	for(i=0;i<comp;i++)
		printf("y%d:%f   gamma%d:%f   fi%d:%f\n",i+1,y[i],i+1,gamma[i],i+1,phi[i]);

break;
}break;


	// DEW T 


	case 4:
		switch(b)
		{
	case 3:


	printf("Please, enter pressure(in kPa):");
	scanf("%lf",&press);

	for(i=0;i<comp-1;i++)
	{
		printf("Please, enter y%d:",i+1);
		scanf("%lf",&y[i]);
		frac=y[i]+frac;
	}
	
	y[i]=1.0-frac;
	for(i=0;i<comp;i++)
		phi[i]=1.0;

	for(i=0;i<comp;i++)
		gamma[i]=1.0;

	for(i=0;i<comp;i++)
		tsat[i]=(-1.0*mtrx[i][8]+(mtrx[i][7]/(-1.0*log(press)+mtrx[i][6])));

	tempc=0.0;
	for(i=0;i<comp;i++)
		tempc+=y[i]*tsat[i];

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));
	
	tempk=tempc+273.15;

	sum1=0.0;
	pjsat=psat[0];

	for(i=0;i<comp;i++)
		sum1+=(y[i]*phi[i]/gamma[i])*(pjsat/psat[i]);

		pjsat=press*sum1;		
		
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{	
			b0=0.083-(0.422/(pow((tempk/tcrit[i][j]),(1.6))));
			b1=0.139-(0.172/(pow((tempk/tcrit[i][j]),(4.2))));
			bij[i][j]=((R*1000.0*tcrit[i][j])/pcrit[i][j])*(b0+wcrit[i][j]*b1);    
		}
	}

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			delta[j][i]=2.0*bij[j][i]-bij[j][j]-bij[i][i];
	}

	sum1=0.0;

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{
			for(k=0;k<comp;k++)
				sum1+=y[j]*y[k]*(2.0*delta[j][i]-delta[j][k]);
		}
		phi[i]=exp((bij[i][i]*(press-psat[i])+0.5*press*sum1)/(R*1000.0*tempk));
		sum1=0.0;
	}
	for(i=0;i<comp;i++)
		x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];

		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);

		sum1=0.0;
	}

	pjsat=psat[0];

	sum1=0.0;
	for(i=0;i<comp;i++)
		sum1+=(y[i]*phi[i]/gamma[i])*(pjsat/psat[i]);
			pjsat=press*sum1;

	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

while(1)
{
	for(i=0;i<comp;i++)
	{
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));
	}
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{	
			b0=0.083-(0.422/(pow((tempk/tcrit[i][j]),(1.6))));
			b1=0.139-(0.172/(pow((tempk/tcrit[i][j]),(4.2))));
			bij[i][j]=((R*1000.0*tcrit[i][j])/pcrit[i][j])*(b0+wcrit[i][j]*b1);    
		}
	}


	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			delta[j][i]=2.0*bij[j][i]-bij[j][j]-bij[i][i];
	}

	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{
			for(k=0;k<comp;k++)
				sum1+=y[j]*y[k]*(2.0*delta[j][i]-delta[j][k]);
		}
		phi[i]=exp((bij[i][i]*(press-psat[i])+0.5*press*sum1)/(R*1000.0*tempk));
		sum1=0.0;
	}

while(1)
{
	for(i=0;i<comp;i++)
		x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

		sum1=0.0;
		for(i=0;i<comp;i++)
			sum1+=x[i];

		for(i=0;i<comp;i++)
			x[i]=x[i]/sum1;
		gammao=(double*)malloc(comp*sizeof(double));
	for(i=0;i<comp;i++)
		gammao[i]=gamma[i];
		
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}
	
	sum1=0.0;
	
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
	
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}
	j=0;
	for(i=0;i<comp;i++)
	{
		if(gammao[i]-gamma[i]<0.1&&gammao[i]-gamma[i]>-0.1)
		j++;
	}
	if(j==comp)
		break;
	}

	pjsat=psat[0];
	sum1=0.0;

	for(i=0;i<comp;i++)
		sum1+=(y[i]*phi[i]/gamma[i])*(pjsat/psat[i]);

	pjsat=press*sum1;
	
	tempco=tempc;
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;
	
	if(tempco-tempc<0.1&&tempco-tempc>-0.1)
	break;
}
	printf("Temperature is %f C\n",tempc);

	for(i=0;i<comp;i++)
		printf("x%d:%f   gamma%d:%f   fi%d:%f\n",i+1,x[i],i+1,gamma[i],i+1,phi[i]);
break;

	case 2:
		printf("Please, enter pressure(in kPa):");
	scanf("%lf",&press);

	for(i=0;i<comp-1;i++)
	{
		printf("Please, enter y%d:",i+1);
		scanf("%lf",&y[i]);
		frac=y[i]+frac;
	}
	
	y[i]=1.0-frac;
	for(i=0;i<comp;i++)
		phi[i]=1.0;

	for(i=0;i<comp;i++)
		gamma[i]=1.0;

	for(i=0;i<comp;i++)
		tsat[i]=(-1.0*mtrx[i][8]+(mtrx[i][7]/(-1.0*log(press)+mtrx[i][6])));

	tempc=0.0;
	for(i=0;i<comp;i++)
		tempc+=y[i]*tsat[i];

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));
	
	tempk=tempc+273.15;

	sum1=0.0;
	pjsat=psat[0];

	for(i=0;i<comp;i++)
		sum1+=(y[i]*phi[i]/gamma[i])*(pjsat/psat[i]);

		pjsat=press*sum1;		
		
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

		for(i=0;i<comp;i++)
		x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];

		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);

		sum1=0.0;
	}

	pjsat=psat[0];

	sum1=0.0;
	for(i=0;i<comp;i++)
		sum1+=(y[i]*phi[i]/gamma[i])*(pjsat/psat[i]);
			pjsat=press*sum1;

	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

while(1)
{
	for(i=0;i<comp;i++)
	{
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));
	}
	

while(1)
{
	for(i=0;i<comp;i++)
		x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

		sum1=0.0;
		for(i=0;i<comp;i++)
			sum1+=x[i];

		for(i=0;i<comp;i++)
			x[i]=x[i]/sum1;
		gammao=(double*)malloc(comp*sizeof(double));
	for(i=0;i<comp;i++)
		gammao[i]=gamma[i];
		
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}
	
	sum1=0.0;
	
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
	
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}
	j=0;
	for(i=0;i<comp;i++)
	{
		if(gammao[i]-gamma[i]<0.1&&gammao[i]-gamma[i]>-0.1)
		j++;
	}
	if(j==comp)
		break;
	}

	pjsat=psat[0];
	sum1=0.0;

	for(i=0;i<comp;i++)
		sum1+=(y[i]*phi[i]/gamma[i])*(pjsat/psat[i]);

	pjsat=press*sum1;
	
	tempco=tempc;
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;
	
	if(tempco-tempc<0.1&&tempco-tempc>-0.1)
	break;
}
	printf("Temperature is %f C\n",tempc);

	for(i=0;i<comp;i++)
		printf("x%d:%f   gamma%d:%f   fi%d:%f\n",i+1,x[i],i+1,gamma[i],i+1,phi[i]);

	break;
	case 1:

		printf("Please, enter pressure(in kPa):");
	scanf("%lf",&press);

	for(i=0;i<comp-1;i++)
	{
		printf("Please, enter y%d:",i+1);
		scanf("%lf",&y[i]);
		frac=y[i]+frac;
	}
	
	y[i]=1.0-frac;
	for(i=0;i<comp;i++)
		phi[i]=1.0;

	for(i=0;i<comp;i++)
		gamma[i]=1.0;

	for(i=0;i<comp;i++)
		tsat[i]=(-1.0*mtrx[i][8]+(mtrx[i][7]/(-1.0*log(press)+mtrx[i][6])));

	tempc=0.0;
	for(i=0;i<comp;i++)
		tempc+=y[i]*tsat[i];

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));
	
	tempk=tempc+273.15;

	sum1=0.0;
	pjsat=psat[0];

	for(i=0;i<comp;i++)
		sum1+=(y[i]*phi[i]/gamma[i])*(pjsat/psat[i]);

		pjsat=press*sum1;		
		
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

		for(i=0;i<comp;i++)
		x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

	
	pjsat=psat[0];

	sum1=0.0;
	for(i=0;i<comp;i++)
		sum1+=(y[i]*phi[i]/gamma[i])*(pjsat/psat[i]);
			pjsat=press*sum1;

	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

while(1)
{
	for(i=0;i<comp;i++)
	{
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));
	}
	

while(1)
{
	for(i=0;i<comp;i++)
		x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

		sum1=0.0;
		for(i=0;i<comp;i++)
			sum1+=x[i];

		for(i=0;i<comp;i++)
			x[i]=x[i]/sum1;
		gammao=(double*)malloc(comp*sizeof(double));
	for(i=0;i<comp;i++)
		gammao[i]=gamma[i];
		
	
	j=0;
	for(i=0;i<comp;i++)
	{
		if(gammao[i]-gamma[i]<0.1&&gammao[i]-gamma[i]>-0.1)
		j++;
	}
	if(j==comp)
		break;
	}

	pjsat=psat[0];
	sum1=0.0;

	for(i=0;i<comp;i++)
		sum1+=(y[i]*phi[i]/gamma[i])*(pjsat/psat[i]);

	pjsat=press*sum1;
	
	tempco=tempc;
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;
	
	if(tempco-tempc<0.1&&tempco-tempc>-0.1)
	break;
}
	printf("Temperature is %f C\n",tempc);

	for(i=0;i<comp;i++)
		printf("x%d:%f   gamma%d:%f   fi%d:%f\n",i+1,x[i],i+1,gamma[i],i+1,phi[i]);

break;
	}
break;

//FLASH 

	case 5:

	switch(b)	{
	case 3:



	printf("Please, enter pressure(in kPa):");
	scanf("%lf",&press);

	pressorg=press;

	z=(double*)malloc(comp * sizeof (double));
	ki=(double*)malloc(comp * sizeof (double));
	
	frac=0.0;

	for(i=0;i<comp-1;i++)
	{
		printf("Please, enter z%d:",i+1);
		scanf("%lf",&z[i]);
		frac=z[i]+frac;
	}
	z[i]=1.0-frac;

	for(i=0;i<comp;i++)
	{
		y[i]=z[i];	
	}

	//dew t 

	y[i]=1.0-frac;
	for(i=0;i<comp;i++)
		phi[i]=1.0;

	for(i=0;i<comp;i++)
		gamma[i]=1.0;

	for(i=0;i<comp;i++)
		tsat[i]=(-1.0*mtrx[i][8]+(mtrx[i][7]/(-1.0*log(press)+mtrx[i][6])));

	tempc=0.0;
	for(i=0;i<comp;i++)
		tempc+=y[i]*tsat[i];

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));
	
	tempk=tempc+273.15;

	sum1=0.0;
	pjsat=psat[0];

	for(i=0;i<comp;i++)
		sum1+=(y[i]*phi[i]/gamma[i])*(pjsat/psat[i]);

		pjsat=press*sum1;		
		
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{	
			b0=0.083-(0.422/(pow((tempk/tcrit[i][j]),(1.6))));
			b1=0.139-(0.172/(pow((tempk/tcrit[i][j]),(4.2))));
			bij[i][j]=((R*1000.0*tcrit[i][j])/pcrit[i][j])*(b0+wcrit[i][j]*b1);    
		}
	}

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			delta[j][i]=2.0*bij[j][i]-bij[j][j]-bij[i][i];
	}

	sum1=0.0;

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{
			for(k=0;k<comp;k++)
				sum1+=y[j]*y[k]*(2.0*delta[j][i]-delta[j][k]);
		}
		phi[i]=exp((bij[i][i]*(press-psat[i])+0.5*press*sum1)/(R*1000.0*tempk));
		sum1=0.0;
	}
	for(i=0;i<comp;i++)
		x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];

		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);

		sum1=0.0;
	}

	pjsat=psat[0];

	sum1=0.0;
	for(i=0;i<comp;i++)
		sum1+=(y[i]*phi[i]/gamma[i])*(pjsat/psat[i]);
			pjsat=press*sum1;

	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

while(1)
{
	for(i=0;i<comp;i++)
	{
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));
	}
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{	
			b0=0.083-(0.422/(pow((tempk/tcrit[i][j]),(1.6))));
			b1=0.139-(0.172/(pow((tempk/tcrit[i][j]),(4.2))));
			bij[i][j]=((R*1000.0*tcrit[i][j])/pcrit[i][j])*(b0+wcrit[i][j]*b1);    
		}
	}


	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			delta[j][i]=2.0*bij[j][i]-bij[j][j]-bij[i][i];
	}

	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{
			for(k=0;k<comp;k++)
				sum1+=y[j]*y[k]*(2.0*delta[j][i]-delta[j][k]);
		}
		phi[i]=exp((bij[i][i]*(press-psat[i])+0.5*press*sum1)/(R*1000.0*tempk));
		sum1=0.0;
	}

while(1)
{
	for(i=0;i<comp;i++)
		x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

		sum1=0.0;
		for(i=0;i<comp;i++)
			sum1+=x[i];

		for(i=0;i<comp;i++)
			x[i]=x[i]/sum1;
		gammao=(double*)malloc(comp*sizeof(double));
	for(i=0;i<comp;i++)
		gammao[i]=gamma[i];
		
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}
	
	sum1=0.0;
	
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
	
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}
	j=0;
	for(i=0;i<comp;i++)
	{
		if(gammao[i]-gamma[i]<0.1&&gammao[i]-gamma[i]>-0.1)
		j++;
	}
	if(j==comp)
		break;
	}

	pjsat=psat[0];
	sum1=0.0;

	for(i=0;i<comp;i++)
		sum1+=(y[i]*phi[i]/gamma[i])*(pjsat/psat[i]);

	pjsat=press*sum1;
	
	tempco=tempc;
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;
	
	if(tempco-tempc<0.1&&tempco-tempc>-0.1)
	break;
}
	tempd=tempc;
	//DEW T 

	//buble t 
	for(i=0;i<comp;i++)
	{
		x[i]=z[i];	
	}
	for(i=0;i<comp;i++)
		phi[i]=1.0;

	for(i=0;i<comp;i++)
		tsat[i]=(-1.0*mtrx[i][8]+(mtrx[i][7]/(-1.0*log(press)+mtrx[i][6])));
	
	tempc=0.0;
	for(i=0;i<comp;i++)
		tempc+=x[i]*tsat[i];
	
	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));
	
	tempk=tempc+273.15;


	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}
	
	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}


	
	pjsat=psat[0];//the first entry of the data file is selected for the calculation of Pjsat
	sum1=0.0;
	pjsato=psat[0];

	for(i=0;i<comp;i++)
		sum1+=(x[i]*gamma[i]/phi[i])*(psat[i]/pjsat);
	
	pjsat=press/sum1;
		
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

while(1)
{
	
	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	for(i=0;i<comp;i++)
		y[i]=x[i]*gamma[i]*psat[i]/(phi[i]*press);

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{	
			b0=0.083-(0.422/(pow((tempk/tcrit[i][j]),(1.6))));
			b1=0.139-(0.172/(pow((tempk/tcrit[i][j]),(4.2))));
			bij[i][j]=((R*1000.0*tcrit[i][j])/pcrit[i][j])*(b0+wcrit[i][j]*b1);    
		}
	}


	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			delta[j][i]=2.0*bij[j][i]-bij[j][j]-bij[i][i];
	}
	sum1=0.0;

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{
			for(k=0;k<comp;k++)
				sum1+=y[j]*y[k]*(2.0*delta[j][i]-delta[j][k]);
		}
		phi[i]=exp((bij[i][i]*(press-psat[i])+0.5*press*sum1)/(R*1000.0*tempk));
		sum1=0.0;
	}

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
	
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}

	for(i=0;i<comp;i++)
	{
		sum1+=(x[i]*gamma[i]/phi[i])*(psat[i]/pjsat);
	}
		pjsat=press/sum1;


	tempco=tempc;
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

	if(tempco-tempc<0.1&&tempco-tempc>-0.1)
	break;
	}
	
	tempb=tempc;
	printf("Tdew is %f C\n",tempd);
	printf("Tbubble is %f C\n",tempb);
	tempc=(tempb+tempd)/2.0;
	printf("Temperature value is %f C\n",tempc);

	tempk=tempc+273.15;

	//bubble t END

	frac=0.0;
	for(i=0;i<comp;i++)
	{
		y[i]=z[i];	
	}

	////dewp 

	for(i=0;i<comp;i++)
		phi[i]=1.0;
	for(i=0;i<comp;i++)
		gamma[i]=1.0;
	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	press=0.0;
	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));

	press=1.0/press;

	for(i=0;i<comp;i++)
			x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
	
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;		
	}
	press=0.0;

	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));
	
	press=1.0/press;

while(1)
{
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{	
			b0=0.083-(0.422/(pow((tempk/tcrit[i][j]),(1.6))));
			b1=0.139-(0.172/(pow((tempk/tcrit[i][j]),(4.2))));
			bij[i][j]=((R*1000.0*tcrit[i][j])/pcrit[i][j])*(b0+wcrit[i][j]*b1);    
		}
	}

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			delta[j][i]=2.0*bij[j][i]-bij[j][j]-bij[i][i];
	}

	sum1=0.0;

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{
			for(k=0;k<comp;k++)
				sum1+=y[j]*y[k]*(2.0*delta[j][i]-delta[j][k]);
		}
		
		phi[i]=exp((bij[i][i]*(press-psat[i])+0.5*press*sum1)/(R*1000.0*tempk));
		sum1=0.0;
	}

	while(1)
	{
		for(i=0;i<comp;i++)
			x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

		gammao=(double*)malloc(comp*sizeof(double));
	for(i=0;i<comp;i++)
		gammao[i]=gamma[i];

		sum1=0.0;
		for(i=0;i<comp;i++)
			sum1+=x[i];

		for(i=0;i<comp;i++)
			x[i]=x[i]/sum1;
		
		for(i=0;i<comp;i++)
		{
			for(j=0;j<comp;j++)
				invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
		}
		sum1=0.0;

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
	
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
	
			sum2+=(x[k]*invv[k][i])/sum3;
		}

		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}
	
	j=0;
	for(i=0;i<comp;i++)
	{
		if(gammao[i]-gamma[i]<0.1&&gammao[i]-gamma[i]>-0.1)
		j++;
	}
	if(j==comp)
		break;
	}
	
	presso=press;
	press=0.0;
	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));
	
	press=1.0/press;
	
	if(press-presso<1.0&&presso-press>-1.0)
	break;
	}
	
	pressdp=press;
	phid=(double*)malloc(comp*sizeof(double));
	gammad=(double*)malloc(comp*sizeof(double));
	for(i=0;i<comp;i++)
	{
		phid[i]=phi[i];
		gammad[i]=gamma[i];
	}




	//bubble p 
	for(i=0;i<comp;i++)
	{
		x[i]=z[i];	
	}
		
	for(i=0;i<comp;i++)
		phi[i]=1.0;

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];

		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
			sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}

	for(i=0;i<comp;i++)
		press+=x[i]*gamma[i]*psat[i]/phi[i];  


while(1)
	{
	for(i=0;i<comp;i++)
		y[i]=x[i]*gamma[i]*psat[i]/(phi[i]*press);

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{	
			b0=0.083-(0.422/(pow((tempk/tcrit[i][j]),(1.6))));
			b1=0.139-(0.172/(pow((tempk/tcrit[i][j]),(4.2))));
			bij[i][j]=((R*1000.0*tcrit[i][j])/pcrit[i][j])*(b0+wcrit[i][j]*b1);
		}
	}

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			delta[j][i]=2.0*bij[j][i]-bij[j][j]-bij[i][i];
	}
	sum1=0.0;

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{
			for(k=0;k<comp;k++)
				sum1+=y[j]*y[k]*(2.0*delta[j][i]-delta[j][k]);
		}
		phi[i]=exp((bij[i][i]*(press-psat[i])+0.5*press*sum1)/(R*1000.0*tempk));
		sum1=0.0;
	}

	presso=press;
	press=0.0;
	for(i=0;i<comp;i++)
		press+=x[i]*gamma[i]*psat[i]/phi[i];  

	if(press-presso<0.1&&presso-press>-0.1)
	break;
}

	pressbp=press;
	phib=(double*)malloc(comp*sizeof(double));
	gammab=(double*)malloc(comp*sizeof(double));

	for(i=0;i<comp;i++)
	{
		phib[i]=phi[i];
		gammab[i]=gamma[i];
	}
//bubble p END
	printf("Pdew is %f kPa at %f K\n",pressdp,tempc);
	printf("Pbubble is %f kPa at %f K\n\n",pressbp,tempc);

	if(pressorg<pressdp||pressorg>pressbp)
	{
		printf("Your pressure value is not between Pdew and Pbuble!\n");
		return -1;
	}

	vi=(pressbp-pressorg)/(pressbp-pressdp);

	est=(pressorg-pressdp)/(pressbp-pressdp);
		
	for(i=0;i<comp;i++)
	{
		phi[i]=(est*(phib[i]-phid[i]))+phid[i];
		gamma[i]=(est*(gammab[i]-gammad[i]))+gammad[i];
	}

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

//    flash 1st loop  
	while(1)
{

	for(i=0;i<comp;i++)
		ki[i]=gamma[i]*psat[i]/(phi[i]*pressorg);
	
	f=0.0;

	for(i=0;i<comp;i++)
		f+=(z[i]*(ki[i]-1.0))/(1.0+vi*(ki[i]-1.0));
	
	fo=f;
	df=0.0;
	for(i=0;i<comp;i++)
		df+=-1.0*(z[i]*(pow((ki[i]-1.0),2.0))/(pow((1.0+vi*(ki[i]-1.0)),2.0)));
	
	vio=vi;
	// calculating V by Newton's Method
	while(1)
	{
	vi=vi-f/df;
	fo=f;
	f=0.0;
	for(i=0;i<comp;i++)
		f+=(z[i]*(ki[i]-1.0))/(1.0+vi*(ki[i]-1.0));

	if(fo-f<0.001&&fo-f>-0.001)
		break;

	df=0.0;
	for(i=0;i<comp;i++)
		df+=-1.0*(z[i]*(pow((ki[i]-1.0),2.0))/(pow((1.0+vi*(ki[i]-1.0)),2.0)));
	}
	
	xo=(double*)malloc(comp*sizeof(double));
	yo=(double*)malloc(comp*sizeof(double));

	for(i=0;i<comp;i++)
	{
		xo[i]=x[i];
		yo[i]=y[i];
		x[i]=z[i]/(1.0+vi*(ki[i]-1.0));
		y[i]=ki[i]*x[i];
	}
	
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	sum1=0.0;
	
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{	
			b0=0.083-(0.422/(pow((tempk/tcrit[i][j]),(1.6))));
			b1=0.139-(0.172/(pow((tempk/tcrit[i][j]),(4.2))));
			bij[i][j]=((R*1000.0*tcrit[i][j])/pcrit[i][j])*(b0+wcrit[i][j]*b1);    
		}
	}

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			delta[j][i]=2.0*bij[j][i]-bij[j][j]-bij[i][i];
	}

	sum1=0.0;

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		{
			for(k=0;k<comp;k++)
				sum1+=y[j]*y[k]*(2.0*delta[j][i]-delta[j][k]);
		}
		phi[i]=exp((bij[i][i]*(press-psat[i])+0.5*press*sum1)/(R*1000.0*tempk));

		sum1=0.0;
	}
	j=0;
	for(i=0;i<comp;i++)
	{
		if(xo[i]-x[i]<0.001&&xo[i]-x[i]>-0.001)
			j++;
		if(yo[i]-y[i]<0.001&&yo[i]-y[i]>-0.001)
			j++;
	}
	if(vio-vi<0.001&&vio-vi>-0.001)
		j++;
	if(j>comp+comp)
		break;
	}
	printf("Vi = %f\n",vi);
	for(i=0;i<comp;i++)
	{
		printf("x%d:%f   y%d:%f gamma%d:%f   fi%d:%f  \n",i+1,x[i],i+1,y[i],i+1,gamma[i],i+1,phi[i]);
	}

break;
	case 2: 

		printf("Please, enter pressure(in kPa):");
	scanf("%lf",&press);

	pressorg=press;

	z=(double*)malloc(comp * sizeof (double));
	ki=(double*)malloc(comp * sizeof (double));
	
	frac=0.0;

	for(i=0;i<comp-1;i++)
	{
		printf("Please, enter z%d:",i+1);
		scanf("%lf",&z[i]);
		frac=z[i]+frac;
	}
	z[i]=1.0-frac;

	for(i=0;i<comp;i++)
	{
		y[i]=z[i];	
	}

	//dew t 

	for(i=0;i<comp;i++)
		phi[i]=1.0;

	for(i=0;i<comp;i++)
		gamma[i]=1.0;

	for(i=0;i<comp;i++)
		tsat[i]=(-1.0*mtrx[i][8]+(mtrx[i][7]/(-1.0*log(press)+mtrx[i][6])));

	tempc=0.0;
	for(i=0;i<comp;i++)
		tempc+=y[i]*tsat[i];

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));
	
	tempk=tempc+273.15;

	sum1=0.0;
	pjsat=psat[0];

	for(i=0;i<comp;i++)
		sum1+=(y[i]*1.0/gamma[i])*(pjsat/psat[i]);

		pjsat=press*sum1;		
		
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

		for(i=0;i<comp;i++)
		x[i]=y[i]*1.0*press/(gamma[i]*psat[i]);

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];

		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);

		sum1=0.0;
	}

	pjsat=psat[0];

	sum1=0.0;
	for(i=0;i<comp;i++)
		sum1+=(y[i]*1.0/gamma[i])*(pjsat/psat[i]);
			pjsat=press*sum1;

	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

while(1)
{
	for(i=0;i<comp;i++)
	{
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));
	}
	

while(1)
{
	for(i=0;i<comp;i++)
		x[i]=y[i]*1.0*press/(gamma[i]*psat[i]);

		sum1=0.0;
		for(i=0;i<comp;i++)
			sum1+=x[i];

		for(i=0;i<comp;i++)
			x[i]=x[i]/sum1;
		gammao=(double*)malloc(comp*sizeof(double));
	for(i=0;i<comp;i++)
		gammao[i]=gamma[i];
		
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}
	
	sum1=0.0;
	
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
	
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}
	j=0;
	for(i=0;i<comp;i++)
	{
		if(gammao[i]-gamma[i]<0.1&&gammao[i]-gamma[i]>-0.1)
		j++;
	}
	if(j==comp)
		break;
	}

	pjsat=psat[0];
	sum1=0.0;

	for(i=0;i<comp;i++)
		sum1+=(y[i]*phi[i]/gamma[i])*(pjsat/psat[i]);

	pjsat=press*sum1;
	
	tempco=tempc;
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;
	
	if(tempco-tempc<0.1&&tempco-tempc>-0.1)
	break;
}
	tempd=tempc;
	//DEW T 

	//buble t 
	for(i=0;i<comp;i++)
		x[i]=z[i];

		for(i=0;i<comp;i++)
		phi[i]=1.0;

	for(i=0;i<comp;i++)
		tsat[i]=(-1.0*mtrx[i][8]+(mtrx[i][7]/(-1.0*log(press)+mtrx[i][6])));
	
	tempc=0.0;
	for(i=0;i<comp;i++)
		tempc+=x[i]*tsat[i];
	
	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));
	
	tempk=tempc+273.15;


	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}
	
	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}


	
	pjsat=psat[0];//the first entry of the data file is selected for the calculation of Pjsat
	sum1=0.0;
	pjsato=psat[0];

	for(i=0;i<comp;i++)
		sum1+=(x[i]*gamma[i]/phi[i])*(psat[i]/pjsat);
	
	pjsat=press/sum1;
		
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

while(1)
{
	
	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	for(i=0;i<comp;i++)
		y[i]=x[i]*gamma[i]*psat[i]/(phi[i]*press);


	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
	
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}

	for(i=0;i<comp;i++)
	{
		sum1+=(x[i]*gamma[i]/phi[i])*(psat[i]/pjsat);
	}
		pjsat=press/sum1;


	tempco=tempc;
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

	if(tempco-tempc<0.1&&tempco-tempc>-0.1)
	break;
	}
	
	tempb=tempc;
	printf("Tdew is %f C\n",tempd);
	printf("Tbubble is %f C\n",tempb);
	tempc=(tempb+tempd)/2.0;
	printf("Temperature value is %f C\n",tempc);

	tempk=tempc+273.15;

	//bubble t 

	frac=0.0;
	for(i=0;i<comp;i++)
	{
		y[i]=z[i];	
	}

	////dewp 

	for(i=0;i<comp;i++)
		phi[i]=1.0;
	for(i=0;i<comp;i++)
		gamma[i]=1.0;
	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	press=0.0;
	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));

	press=1.0/press;

	for(i=0;i<comp;i++)
			x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
	
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;		
	}
	press=0.0;

	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));
	
	press=1.0/press;

while(1)
{
	while(1)
	{
		for(i=0;i<comp;i++)
			x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

		gammao=(double*)malloc(comp*sizeof(double));
	for(i=0;i<comp;i++)
		gammao[i]=gamma[i];

		sum1=0.0;
		for(i=0;i<comp;i++)
			sum1+=x[i];

		for(i=0;i<comp;i++)
			x[i]=x[i]/sum1;
		
		for(i=0;i<comp;i++)
		{
			for(j=0;j<comp;j++)
				invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
		}
		sum1=0.0;

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
	
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
	
			sum2+=(x[k]*invv[k][i])/sum3;
		}

		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}
	
	j=0;
	for(i=0;i<comp;i++)
	{
		if(gammao[i]-gamma[i]<0.1&&gammao[i]-gamma[i]>-0.1)
		j++;
	}
	if(j==comp)
		break;
	}
	
	presso=press;
	press=0.0;
	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));
	
	press=1.0/press;
	
	if(press-presso<1.0&&presso-press>-1.0)
	break;
	}
	
	pressdp=press;
	phid=(double*)malloc(comp*sizeof(double));
	gammad=(double*)malloc(comp*sizeof(double));
	for(i=0;i<comp;i++)
	{
		phid[i]=phi[i];
		gammad[i]=gamma[i];
	}




	//bubble p 

	for(i=0;i<comp;i++)
		x[i]=z[i];
	for(i=0;i<comp;i++)
		phi[i]=1.0;

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
		invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];

		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
			sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}

	for(i=0;i<comp;i++)
		press+=x[i]*gamma[i]*psat[i]/phi[i];  


while(1)
	{
	for(i=0;i<comp;i++)
		y[i]=x[i]*gamma[i]*psat[i]/(phi[i]*press);

	presso=press;
	press=0.0;
	for(i=0;i<comp;i++)
		press+=x[i]*gamma[i]*psat[i]/phi[i];  

	if(press-presso<0.1&&presso-press>-0.1)
	break;
	}

	pressbp=press;
	phib=(double*)malloc(comp*sizeof(double));
	gammab=(double*)malloc(comp*sizeof(double));

	for(i=0;i<comp;i++)
	{
		phib[i]=1.0;
		gammab[i]=gamma[i];
	}
//bubble p 
	printf("Pdew is %f kPa at %f K\n",pressdp,tempc);
	printf("Pbubble is %f kPa at %f K\n\n",pressbp,tempc);

	if(pressorg<pressdp||pressorg>pressbp)
	{
		printf("Your pressure value is not between Pdew and Pbuble!\n");
		return -1;
	}

	vi=(pressbp-pressorg)/(pressbp-pressdp);

	est=(pressorg-pressdp)/(pressbp-pressdp);
		
	for(i=0;i<comp;i++)
	{
		phi[i]=1.0;
		gamma[i]=(est*(gammab[i]-gammad[i]))+gammad[i];
	}

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

//    flash 1st loop  
	while(1)
{

	for(i=0;i<comp;i++)
		ki[i]=gamma[i]*psat[i]/(phi[i]*pressorg);
	
	f=0.0;

	for(i=0;i<comp;i++)
		f+=(z[i]*(ki[i]-1.0))/(1.0+vi*(ki[i]-1.0));
	
	fo=f;
	df=0.0;
	for(i=0;i<comp;i++)
		df+=-1.0*(z[i]*(pow((ki[i]-1.0),2.0))/(pow((1.0+vi*(ki[i]-1.0)),2.0)));
	
	vio=vi;
	// calculating V by Newton's Method
	while(1)
	{
	vi=vi-f/df;
	fo=f;
	f=0.0;
	for(i=0;i<comp;i++)
		f+=(z[i]*(ki[i]-1.0))/(1.0+vi*(ki[i]-1.0));

	if(fo-f<0.001&&fo-f>-0.001)
		break;

	df=0.0;
	for(i=0;i<comp;i++)
		df+=-1.0*(z[i]*(pow((ki[i]-1.0),2.0))/(pow((1.0+vi*(ki[i]-1.0)),2.0)));
	}
	
	xo=(double*)malloc(comp*sizeof(double));
	yo=(double*)malloc(comp*sizeof(double));

	for(i=0;i<comp;i++)
	{
		xo[i]=x[i];
		yo[i]=y[i];
		x[i]=z[i]/(1.0+vi*(ki[i]-1.0));
		y[i]=ki[i]*x[i];
	}
	
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			invv[i][j]=(mtrx[j][5]/mtrx[i][5])*exp((-1.0*mtrx[i][j+9])/(R*tempk));
	}

	sum1=0.0;
	
	for(i=0;i<comp;i++)
	{
		for(j=0;j<comp;j++)
			sum1+=x[j]*invv[i][j];
		sum2=0.0;
		for(k=0;k<comp;k++)
		{
		    sum3=0.0;
			for(j=0;j<comp;j++)
				sum3+=x[j]*invv[k][j];
			sum2+=(x[k]*invv[k][i])/sum3;
		}
		gamma[i]=exp(1.0-log(sum1)-sum2);
		sum1=0.0;
	}

	
	j=0;
	for(i=0;i<comp;i++)
	{
		if(xo[i]-x[i]<0.001&&xo[i]-x[i]>-0.001)
			j++;
		if(yo[i]-y[i]<0.001&&yo[i]-y[i]>-0.001)
			j++;
	}
	if(vio-vi<0.001&&vio-vi>-0.001)
		j++;
	if(j>comp+comp)
		break;
	}
	printf("Vi = %f\n",vi);
	for(i=0;i<comp;i++)
	{
		printf("x%d:%f   y%d:%f gamma%d:%f   fi%d:%f  \n",i+1,x[i],i+1,y[i],i+1,gamma[i],i+1,phi[i]);
	}

break;
	case 1:      /////////////////// ////////////////////////////////////////////////////////////  case 11cse11111case11111111111111case111111

		printf("Please, enter pressure(in kPa):");
	scanf("%lf",&press);

	pressorg=press;

	z=(double*)malloc(comp * sizeof (double));
	ki=(double*)malloc(comp * sizeof (double));
	
	frac=0.0;

	for(i=0;i<comp-1;i++)
	{
		printf("Please, enter z%d:",i+1);
		scanf("%lf",&z[i]);
		frac=z[i]+frac;
	}
	z[i]=1.0-frac;

	for(i=0;i<comp;i++)
		y[i]=z[i];

	for(i=0;i<comp;i++)
		phi[i]=1.0;

	for(i=0;i<comp;i++)
		gamma[i]=1.0;

	for(i=0;i<comp;i++)
		tsat[i]=(-1.0*mtrx[i][8]+(mtrx[i][7]/(-1.0*log(press)+mtrx[i][6])));

	tempc=0.0;
	for(i=0;i<comp;i++)
		tempc+=y[i]*tsat[i];

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));
	
	tempk=tempc+273.15;

	sum1=0.0;
	pjsat=psat[0];

	for(i=0;i<comp;i++)
		sum1+=(y[i]*phi[i]/gamma[i])*(pjsat/psat[i]);

		pjsat=press*sum1;		
		
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

		for(i=0;i<comp;i++)
		x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

	
	pjsat=psat[0];

	sum1=0.0;
	for(i=0;i<comp;i++)
		sum1+=(y[i]*phi[i]/gamma[i])*(pjsat/psat[i]);
			pjsat=press*sum1;

	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;

while(1)
{
	for(i=0;i<comp;i++)
	{
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));
	}
	

while(1)
{
	for(i=0;i<comp;i++)
		x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

		sum1=0.0;
		for(i=0;i<comp;i++)
			sum1+=x[i];

		for(i=0;i<comp;i++)
			x[i]=x[i]/sum1;
		gammao=(double*)malloc(comp*sizeof(double));
	for(i=0;i<comp;i++)
		gammao[i]=gamma[i];
		
	
	j=0;
	for(i=0;i<comp;i++)
	{
		if(gammao[i]-gamma[i]<0.1&&gammao[i]-gamma[i]>-0.1)
		j++;
	}
	if(j==comp)
		break;
	}

	pjsat=psat[0];
	sum1=0.0;

	for(i=0;i<comp;i++)
		sum1+=(y[i]*phi[i]/gamma[i])*(pjsat/psat[i]);

	pjsat=press*sum1;
	
	tempco=tempc;
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;
	
	if(tempco-tempc<0.1&&tempco-tempc>-0.1)
	break;
}	tempd=tempc;
	//DEW T 

	//buble t 

for(i=0;i<comp;i++)
x[i]=z[i];
	for(i=0;i<comp;i++)
	{
		phi[i]=1.0;
		gamma[i]=1.0;
	}

	
	for(i=0;i<comp;i++)
		tsat[i]=(-1.0*mtrx[i][8]+(mtrx[i][7]/(-1.0*log(press)+mtrx[i][6])));	
	tempc=0.0;
	for(i=0;i<comp;i++)
		tempc+=x[i]*tsat[i];

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	tempk=tempc+273.15;


	pjsat=psat[0];//the first entry of the data file is selected for the calculation of Pjsat
	sum1=0.0;
	pjsato=psat[0];

	for(i=0;i<comp;i++)
		sum1+=(x[i]*gamma[i]/phi[i])*(psat[i]/pjsat);
	
	pjsat=press/sum1;
		
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;


while(1)
{
	
	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));


	for(i=0;i<comp;i++)
		y[i]=x[i]*gamma[i]*psat[i]/(phi[i]*press);




	sum1=0.0;
	for(i=0;i<comp;i++)
	{
		sum1+=(x[i]*gamma[i]/phi[i])*(psat[i]/pjsat);
	}
		pjsat=press/sum1;


	tempco=tempc;
	tempc=(-1.0*mtrx[0][8]+(mtrx[0][7]/(-1.0*log(pjsat)+mtrx[0][6])));
	tempk=tempc+273.15;


	if(tempco-tempc<0.1&&tempco-tempc>-0.1)
	break;
	}
	
	tempb=tempc;
	printf("Tdew is %f C\n",tempd);
	printf("Tbubble is %f C\n",tempb);
	tempc=(tempb+tempd)/2.0;
	printf("Temperature value is %f C\n",tempc);

	tempk=tempc+273.15;

	//bubble t 

	frac=0.0;
	for(i=0;i<comp;i++)
	{
		y[i]=z[i];	
	}

	////dewp 

	for(i=0;i<comp;i++)
		phi[i]=1.0;
	for(i=0;i<comp;i++)
		gamma[i]=1.0;
	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	press=0.0;
	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));

	press=1.0/press;

	for(i=0;i<comp;i++)
			x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

	press=0.0;

	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));
	
	press=1.0/press;

while(1)
{
	while(1)
	{
		for(i=0;i<comp;i++)
			x[i]=y[i]*phi[i]*press/(gamma[i]*psat[i]);

		gammao=(double*)malloc(comp*sizeof(double));
	for(i=0;i<comp;i++)
		gammao[i]=gamma[i];

		sum1=0.0;
		for(i=0;i<comp;i++)
			sum1+=x[i];

		for(i=0;i<comp;i++)
			x[i]=x[i]/sum1;
			
	j=0;
	for(i=0;i<comp;i++)
	{
		if(gammao[i]-gamma[i]<0.1&&gammao[i]-gamma[i]>-0.1)
		j++;
	}
	if(j==comp)
		break;
	}
	
	presso=press;
	press=0.0;
	for(i=0;i<comp;i++)
		press+=((y[i]*phi[i])/(gamma[i]*psat[i]));
	
	press=1.0/press;
	
	if(press-presso<1.0&&presso-press>-1.0)
	break;
	}
	
	pressdp=press;
	phid=(double*)malloc(comp*sizeof(double));
	gammad=(double*)malloc(comp*sizeof(double));
	for(i=0;i<comp;i++)
	{
		phid[i]=1.0;
		gammad[i]=1.0;
		gamma[i]=1.0;
		phi[i]=1.0;
	}




	//bubble p 
for(i=0;i<comp;i++)
		x[i]=z[i];
	

	for(i=0;i<comp;i++)
		phi[i]=1.0;
	for(i=0;i<comp;i++)
		gamma[i]=1.0;

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

	for(i=0;i<comp;i++)
		press+=x[i]*gamma[i]*psat[i]/phi[i];  


while(1)
	{
	for(i=0;i<comp;i++)
		y[i]=x[i]*gamma[i]*psat[i]/(phi[i]*press);

	presso=press;
	press=0.0;
	for(i=0;i<comp;i++)
		press+=x[i]*gamma[i]*psat[i]/phi[i];  

	if(press-presso<0.1&&presso-press>-0.1)
	break;
	}

	pressbp=press;
	phib=(double*)malloc(comp*sizeof(double));
	gammab=(double*)malloc(comp*sizeof(double));

	for(i=0;i<comp;i++)
	{
		phib[i]=1.0;
		gammab[i]=1.0;
		gamma[i]=1.0;
	}
//bubble p 
	printf("Pdew is %f kPa at %f C\n",pressdp,tempc);
	printf("Pbubble is %f kPa at %f C,\n\n",pressbp,tempc);

	if(pressorg<pressdp||pressorg>pressbp)
	{
		printf("Your pressure value is not between Pdew and Pbuble!\n");
		return -1;
	}

	vi=(pressbp-pressorg)/(pressbp-pressdp);

	est=(pressorg-pressdp)/(pressbp-pressdp);
		
	for(i=0;i<comp;i++)
	{
		phi[i]=1.0;
		gamma[i]=1.0;
	}

	for(i=0;i<comp;i++)
		psat[i]=exp(mtrx[i][6]-(mtrx[i][7]/(tempc+mtrx[i][8])));

//    flash 1st loop  
	while(1)
{

	for(i=0;i<comp;i++)
		ki[i]=gamma[i]*psat[i]/(phi[i]*pressorg);
	
	f=0.0;

	for(i=0;i<comp;i++)
		f+=(z[i]*(ki[i]-1.0))/(1.0+vi*(ki[i]-1.0));
	
	fo=f;
	df=0.0;
	for(i=0;i<comp;i++)
		df+=-1.0*(z[i]*(pow((ki[i]-1.0),2.0))/(pow((1.0+vi*(ki[i]-1.0)),2.0)));
	
	vio=vi;
	// calculating V by Newton's Method
	while(1)
	{
	vi=vi-f/df;
	fo=f;
	f=0.0;
	for(i=0;i<comp;i++)
		f+=(z[i]*(ki[i]-1.0))/(1.0+vi*(ki[i]-1.0));

	if(fo-f<0.001&&fo-f>-0.001)
		break;

	df=0.0;
	for(i=0;i<comp;i++)
		df+=-1.0*(z[i]*(pow((ki[i]-1.0),2.0))/(pow((1.0+vi*(ki[i]-1.0)),2.0)));
	}
	
	xo=(double*)malloc(comp*sizeof(double));
	yo=(double*)malloc(comp*sizeof(double));

	for(i=0;i<comp;i++)
	{
		xo[i]=x[i];
		yo[i]=y[i];
		x[i]=z[i]/(1.0+vi*(ki[i]-1.0));
		y[i]=ki[i]*x[i];
	}
	


	sum1=0.0;
	

	
	j=0;
	for(i=0;i<comp;i++)
	{
		if(xo[i]-x[i]<0.001&&xo[i]-x[i]>-0.001)
			j++;
		if(yo[i]-y[i]<0.001&&yo[i]-y[i]>-0.001)
			j++;
	}
	if(vio-vi<0.001&&vio-vi>-0.001)
		j++;
	if(j>comp+comp)
		break;
	}
	printf("Vi = %f\n",vi);
	for(i=0;i<comp;i++)
	{
		printf("x%d:%f   y%d:%f gamma%d:%f   fi%d:%f  \n",i+1,x[i],i+1,y[i],i+1,gamma[i],i+1,phi[i]);
	}

break;	}
	break;
	}


getch();
	fclose(input);
	return 7;
}