//======================================================
// java Termodinamik paketi
// Class sogut soğutma akışkanların özellikleri
// Dr. Turhan Coban
// Ege Üniversitesi Mühendislik fakultesi
// Makina Bölümü
// 
//======================================================
// Dosya adı   : sogut.java
// Bu dosya soğutucu akışkanların termodinamik özelliklerini 
// Martin-Hou hal denklemiyle hesaplar
// 
// =====================================================
// P : pressure as a function of temperature and specific volume
//     equation of state
// Pi : pressure as a function of temperature and specific volume
//      also predict saturation presure if the volume in saturation region
// dPdv : derivative of P with respect to v
// dPdt : derivative of P with respect to t
// property : method to calculate several thermodynamic properties
// u  : internal energy as a function of t and v
// h  : enthalpy as a function of t and v
// s  : entropy as a function of t and v
// Ps    : saturated vapor pressure as a function of T
// Ts    : saturated vapor temperature as a function of T
// vsl_t :saturated liquid specific volume as a function of T
// vsg_t :saturated vapor specific volume as a function of T
// v_tp  : specific volume as a function of T and v
// v_th  : specific volume as a function of T and h
// v_tu  : specific volume as a function of T and u
// v_tx  : specific volume as a function of T and x
// v_ts  : specific volume as a function of T and s
// t_pv  : temperature as a function of pressure and specific volume
// t_ph  : temperature as a function of pressure and enthalpy
// t_pu  : temperature as a function of pressure and internal energy
// t_ps  : temperature as a function of pressure and internal entropy
// x_tv : quality as a function of temperature and specific volume
// toString : result output through a string variable

import java.io.*;
import java.net.*;
import java.util.*;

class sogut1
{
	
boolean SI=true;   //switch for SI/EN units
public sogut1Data rd;
final String u[]={"SI","EN"};
String unit=u[0];
//=======constroctors=======
public sogut1()
{ readsogut("R22");}

public sogut1(String rName)
{ readsogut(rName);}

public String readGasNames()
{ String s1=rd.refrigerantName;String s2=readsogut("null");readsogut(s1);return s2;}
//property calculation methods===========
public double[] property_EN(String s, double v1, double v2)
{
if(s.charAt(0)=='t') v1=(v1-32.0)/1.8;
else if(s.charAt(0)=='p') v1/=145.04;
else if(s.charAt(0)=='v') v1/=16.01877892;
if(s.charAt(1)=='t') v2=(v2-32)/1.8;
else if(s.charAt(1)=='p') v2/=14.504;
else if(s.charAt(1)=='v') v2/=16.01877892;
else if(s.charAt(1)=='h') v2/=0.42992;
else if(s.charAt(1)=='u') v2/=0.42992;
else if(s.charAt(1)=='s') v2/=0.2388444444;
double pp[]=property_SI(s,v1,v2);
pp[0]*=145.04;
pp[1]=pp[1]*1.8+32;
pp[2]*=16.01877892;
pp[3]*=0.42992;
pp[4]*=0.42992;
pp[5]*=0.2388444444;
pp[7]=1.0/pp[2];
return pp;
}

public double[] property_SI(String s, double v1, double v2)
{
//---- cases of String s --------
//    tv tp th tu ts tx
// pv    pt ph pu ps px
// vp    vt

//---- Exit properties ---------
// r[0] P Pressure MPa
// r[1] t temperature degree C
// r[2] v specific volume m^3/kg
// r[3] h KJ/kg 
// r[4] u KJ/kg
// r[5] s KJ/kgK
// r[6] x kg vapor/kg total phase
// r[7] ro density kg/m^3

double[] r=new double[8];
char s1=s.charAt(0);
char s2=s.charAt(1);
r[6]=x_tv(v1,v2);
if(s.equals("tv")) //tv case
{
 r[6]=x_tv(v1,v2);
 if(r[6]>=0 && r[6]<=1)
 { r[0]=Ps(v1);
   double vsl=vsl_t(v1);
   double vsg=vsg_t(v1);
   r[3]=h(v1,vsl)*(1-r[6])+ h(v1,vsl)*r[6];
   r[4]=u(v1,vsl)*(1-r[6])+ u(v1,vsl)*r[6];
   r[5]=s(v1,vsl)*(1-r[6])+ s(v1,vsl)*r[6];
 }
 else
 { r[0]=P(v1,v2);
   r[3]=h(v1,v2);
   r[4]=u(v1,v2);
   r[5]=s(v1,v2);
 }
   r[1]=v1;
   r[2]=v2;
} //end of tv case
else if(s.equals("pv"))  // pv case
{
 r[6]=x_pv(v1,v2);
 if(r[6]>=0 && r[6]<=1)
 { r[1]=Ts(v1);
   double vsl=vsl_t(r[1]);
   double vsg=vsg_t(r[1]);
   r[3]=h(r[1],vsl)*(1-r[6])+ h(r[1],vsl)*r[6];
   r[4]=u(r[1],vsl)*(1-r[6])+ u(r[1],vsl)*r[6];
   r[5]=s(r[1],vsl)*(1-r[6])+ s(r[1],vsl)*r[6];
 }
 else
 {
   r[1]=t_pv(v1,v2);
   r[3]=h(r[1],v2);
   r[4]=u(r[1],v2);
   r[5]=s(r[1],v2);
 }
  r[0]=v1;
  r[2]=v2;
} //end of pv case
else if(s.equals("tp")) //tp case
{
   r[0]=v2;
   r[1]=v1;
   r[2]=v_tp(v1,v2);
   r[3]=h(r[1],r[2]);
   r[4]=u(r[1],r[2]);
   r[5]=s(r[1],r[2]);
   r[6]=x_tv(r[1],r[2]);
} //end of tp case
else if(s.equals("th"))
{
 double hsl=h(v1,vsl_t(v1));
 double hsg=h(v1,vsg_t(v1));
 double vth=v_th(v1,v2);
   r[0]=P(v1,vth);
   r[1]=v1;
   r[2]=vth;
   r[3]=v2;
   r[4]=u(v1,vth);
   r[5]=s(v1,vth);
   r[6]=x_tv(v1,r[2]);
 if(r[6]>=0 && r[6]<=1)
 { r[0]=Ps(v1);
   r[1]=v1;
   double vsl=vsl_t(v1);
   double vsg=vsg_t(v1);
   r[2]=vsl*(1-r[6])+ vsg*r[6];
   r[3]=v2;
   r[4]=h(v1,vsl)*(1-r[6])+ h(v1,vsg)*r[6];
   r[5]=s(v1,vsl)*(1-r[6])+ s(v1,vsg)*r[6];
  }
} //end of th case
else if(s.equals("tu"))
{
   double vtu=v_tu(v1,v2);
   r[0]=P(v1,vtu);
   r[1]=v1;
   r[2]=vtu;
   r[3]=h(v1,vtu);
   r[4]=v2;
   r[5]=s(v1,vtu);
   r[6]=x_tv(v1,r[2]);
   if(r[6]>=0 && r[6]<=1)
   {
      r[0]=Ps(v1);
      r[1]=v1;
      double vsl=vsl_t(v1);
      double vsg=vsg_t(v1);
      r[2]=vsl*(1-r[6])+ vsg*r[6];
      r[3]=h(v1,vsl)*(1-r[6])+ h(v1,vsg)*r[6];
      r[4]=v2;
      r[5]=s(v1,vsl)*(1-r[6])+ s(v1,vsg)*r[6];
  }
} // end of tu case
else if(s.equals("ts"))
{
   double vts=v_ts(v1,v2);
   r[0]=P(v1,vts);
   r[1]=v1;
   r[2]=v_ts(v1,v2);
   r[3]=h(v1,vts);
   r[4]=u(v1,vts);
   r[5]=v2;
   r[6]=x_tv(v1,r[2]);
 if(r[6]>=0 && r[6]<=1)
 { r[0]=Ps(v1);
   double vsl=vsl_t(v1);
   double vsg=vsg_t(v1);
   r[2]=vsl*(1-r[6])+ vsg*r[6];
   r[3]=h(v1,vsl)*(1-r[6])+ h(v1,vsg)*r[6];
   r[4]=u(v1,vsl)*(1-r[6])+ u(v1,vsg)*r[6];
 }

} //end of ts case
else if(s.equals("tx"))
{
   r[6]=v2;
   r[0]=Ps(v1);
   r[1]=v1;
   double hl=hsl_t(v1);
   double hg=hsg_t(v1);
   double vsl=vsl_t(v1);
   double vsg=vsg_t(v1);
   double ul=hl-r[0]*vsl;
   double ug=hg-r[0]*vsg;
   r[2]=vsl*(1-r[6])+ vsg*r[6];
   r[3]=hl*(1-r[6])+ hg*r[6];
   r[4]=ul*(1-r[6])+ ug*r[6];
   r[5]=ssl_t(v1)*(1-r[6])+ ssg_t(v1)*r[6];
} //end of tx case
else if(s.equals("vt"))
{
 r=property_SI("tv",v2,v1);
} //end of vt case
else if(s.equals("vp"))
{
 r=property_SI("pv",v2,v1);
} //end of vp case
else if(s.equals("pt"))
{
 r=property_SI("tp",v2,v1);
} //end of tp case
else if(s.equals("ph"))
{
 double ts=Ts(v1);
 double vsl=vsl_t(ts);
 double vsg=vsg_t(ts);
 double hsl=h(ts,vsl);
 double hsg=h(ts,vsg);
 double vph=v_ph(v1,v2);
   r[0]=v1;
   r[1]=t_ph(v1,v2);
   r[2]=vph;
   r[3]=v2;
   r[4]=u(r[1],vph);
   r[5]=s(r[1],vph);
   if( v1 < rd.Pc)
   { r[6]=(v2-hsl)/(hsg-hsl);}
   else r[6]=3;
 if(r[6]>=0 && r[6]<=1)
 { r[0]=v1;
   r[1]=ts;
   r[2]=vsl*(1-r[6])+ vsg*r[6];
   r[3]=v2;
   r[4]=u(ts,vsl)*(1-r[6])+ u(ts,vsg)*r[6];
   r[5]=s(ts,vsl)*(1-r[6])+ s(ts,vsg)*r[6];
  }
} //end of ph case
else if(s.equals("pu"))
{
 double ts=Ts(v1);
 double vsl=vsl_t(ts);
 double vsg=vsg_t(ts);
 double usl=u(ts,vsl);
 double usg=u(ts,vsg);
 double vpu=v_pu(v1,v2);
   r[0]=v1;
   r[1]=t_pu(v1,v2);
   r[2]=vpu;
   r[3]=h(r[1],r[2]);
   r[4]=v2;
   r[5]=s(r[1],r[2]);
   if( v1 < rd.Pc)
   { r[6]=(v2-usl)/(usg-usl);}
   else r[6]=3;
 if(r[6]>=0 && r[6]<=1)
 { r[0]=v1;
   r[1]=ts;
   r[2]=vsl*(1-r[6])+ vsg*r[6];
   r[3]=h(ts,vsl)*(1-r[6])+ h(ts,vsg)*r[6];
   r[4]=v2;
   r[5]=s(ts,vsl)*(1-r[6])+ s(ts,vsg)*r[6];
  }
} //end of pu case
else if(s.equals("ps"))
{
 double ts=Ts(v1);
 double vsl=vsl_t(ts);
 double vsg=vsg_t(ts);
 double ssl=s(ts,vsl);
 double ssg=s(ts,vsg);
 double vps=v_ps(v1,v2);
   r[0]=v1;
   r[1]=t_ps(v1,v2);
   r[2]=vps;
   r[3]=h(r[1],r[2]);
   r[4]=u(r[1],r[2]);
   r[5]=v2;
   if( v1 < rd.Pc)
   { r[6]=(v2-ssl)/(ssg-ssl);}
   else r[6]=3;
 if(r[6]>=0 && r[6]<=1)
 { r[0]=v1;
   r[1]=ts;
   r[2]=vsl*(1-r[6])+ vsg*r[6];
   r[3]=h(ts,vsl)*(1-r[6])+ h(ts,vsg)*r[6];
   r[4]=u(ts,vsl)*(1-r[6])+ u(ts,vsg)*r[6];
   r[5]=v2;
  }
} //end of ps case
else if(s.equals("px"))
{
   r[6]=v2;
   double ts=Ts(v1);
   r[0]=v1;
   r[1]=ts;
   double hl=hsl_t(ts);
   double hg=hsg_t(ts);
   double vsl=vsl_t(ts);
   double vsg=vsg_t(ts);
   double ul=hl-r[0]*vsl;
   double ug=hg-r[0]*vsg;
   r[2]=vsl*(1-r[6])+ vsg*r[6];
   r[3]=hl*(1-r[6])+ hg*r[6];
   r[4]=ul*(1-r[6])+ ug*r[6];
   r[5]=ssl_t(v1)*(1-r[6])+ ssg_t(v1)*r[6];
} //end of px case
r[7]=1.0/r[2];
return r;
}

public double[] property(String s, double v1, double v2)
{
if(unit.equals(u[0])) return property_SI(s,v1,v2);
else                 return property_EN(s,v1,v2);
}
// ===== Equation of State P(t,v) =======
public double P(double t, double v)
{
//basic equation of state	
//v m^3/kg
double T=t+273.15;
double Tr=T/rd.Tc;
double V=v*rd.M;
double ror=rd.Vc/V;
double z=0;
double Tr2=Tr*Tr;double Tr3=Tr2*Tr;double Tr4=Tr2*Tr2;double Tr6=Tr2*Tr4;
double ror2=ror*ror;double ror3=ror2*ror;double ror4=ror2*ror2;double ror5=ror3*ror2;
double ror6=ror4*ror2;double ror7=ror6*ror;double ror8=ror6*ror2;double ror9=ror8*ror;
double ror10=ror8*ror2;double ror11=ror10*ror;
double B[]=new double[11];
double R=rd.R/rd.M;
if(rd.refrigerantName.equals("R22"))  
{B[0]=(rd.A[0]+rd.A[1]/Tr+rd.A[2]/Tr4+rd.A[3]/Tr6);
 B[1]=(rd.A[4]+rd.A[5]/Tr+rd.A[6]/Tr4);
 B[2]=(rd.A[7]+rd.A[8]/Tr+rd.A[9]/Tr2+rd.A[10]/Tr4);
 B[3]=(rd.A[11]/Tr+rd.A[12]/Tr4);
 B[4]=(rd.A[13]/Tr2+rd.A[14]/Tr3+rd.A[15]/Tr4+rd.A[16]/Tr6);
 B[5]=(rd.A[17]);
 B[6]=(rd.A[18]/Tr2+rd.A[19]/Tr6);
 B[7]=(rd.A[20]/Tr);
 B[8]=(rd.A[21]/Tr2+rd.A[22]/Tr6);
 B[9]=(rd.A[23]/Tr+rd.A[24]/Tr6);
 B[10]=(rd.A[25]/Tr+rd.A[26]/Tr2+rd.A[27]/Tr6);
}
else{
 B[0]=(rd.A[0]+rd.A[1]/Tr+rd.A[2]/Tr2+rd.A[3]/Tr3);
 B[1]=(rd.A[4]+rd.A[5]/Tr+rd.A[6]/Tr2);
 B[2]=(rd.A[7]+rd.A[8]/Tr+rd.A[9]/Tr2);
 B[3]=(rd.A[10]+rd.A[11]/Tr);
 B[4]=0.0;
 B[5]=(rd.A[12]+rd.A[13]/Tr);
 B[6]=0.0;
 B[7]=(rd.A[14]);
 B[8]=0.0;
 B[9]=0.0;
 B[10]=0.0;	
}
z=1.0+B[0]*ror+B[1]*ror2+B[2]*ror3+B[3]*ror4+B[4]*ror5+B[5]*ror6+B[6]*ror7+B[7]*ror8+
           B[8]*ror9+B[9]*ror10+B[10]*ror11; 
return z*R*T/v;
}

double Pi(double t, double v)
{
if(vsl_t(t)<=v && v>=vsg_t(t)) {return Ps(t);}
return P(t,v);
}

//==== derivatives of equation of state ====
public double dPdT(double t, double v)
{
double T=t+273.15;
double Tr=T/rd.Tc;
double V=v*rd.M;
double ror=rd.Vc/V;
double dzdTr=0;
double Tr2=Tr*Tr;double Tr3=Tr2*Tr;double Tr4=Tr2*Tr2;double Tr5=Tr4*Tr;
double Tr6=Tr2*Tr4;double Tr7=Tr6*Tr;
double ror2=ror*ror;double ror3=ror2*ror;double ror4=ror2*ror2;double ror5=ror3*ror2;
double ror6=ror4*ror2;double ror7=ror6*ror;double ror8=ror6*ror2;double ror9=ror8*ror;
double ror10=ror8*ror2;double ror11=ror10*ror;
double dB[]=new double[11];
double R=rd.R/rd.M;
if(rd.refrigerantName.equals("R22"))  
{dB[0]=(-rd.A[1]/Tr2-4.0*rd.A[2]/Tr5-6.0*rd.A[3]/Tr7);
 dB[1]=(-rd.A[5]/Tr2-4.0*rd.A[6]/Tr5);
 dB[2]=(-rd.A[8]/Tr2-2.0*rd.A[9]/Tr3-4.0*rd.A[10]/Tr5);
 dB[3]=(-rd.A[11]/Tr2-4.0*rd.A[12]/Tr5);
 dB[4]=(-2.0*rd.A[13]/Tr3-3.0*rd.A[14]/Tr4-4.0*rd.A[15]/Tr5-6.0*rd.A[16]/Tr7);
 dB[5]=0.0;
 dB[6]=(-2.0*rd.A[18]/Tr3-6.0*rd.A[19]/Tr7);
 dB[7]=(-rd.A[20]/Tr2);
 dB[8]=(-2.0*rd.A[21]/Tr3-6.0*rd.A[22]/Tr7);
 dB[9]=(-rd.A[23]/Tr2-6.0*rd.A[24]/Tr7);
 dB[10]=(-rd.A[25]/Tr2-2.0*rd.A[26]/Tr3-6.0*rd.A[27]/Tr7);	
}
else{
 dB[0]=(-rd.A[1]/Tr2-2.0*rd.A[2]/Tr3-3.0*rd.A[3]/Tr4);
 dB[1]=(-rd.A[5]/Tr2-2.0*rd.A[6]/Tr3);
 dB[2]=(-rd.A[8]/Tr2-2.0*rd.A[9]/Tr3);
 dB[3]=(-rd.A[11]/Tr2);
 dB[4]=0.0;
 dB[5]=(-2.0*rd.A[13]/Tr2);
 dB[6]=0.0;
 dB[7]=0.0;
 dB[8]=0.0;
 dB[9]=0.0;
 dB[10]=0.0;
}
dzdTr= dB[0]*ror+dB[1]*ror2+dB[2]*ror3+dB[3]*ror4+dB[4]*ror5+dB[6]*ror7+dB[7]*ror8+dB[8]*ror9+
    dB[9]*ror10+dB[10]*ror11; 
return dzdTr*R*Tr/v+P(t,v)/T;
}

public  double dPdT1(double x,double v)
{
// numerically calculated dP/dT
double h0=0.0256;
int i,m;
int n=7;
//derivative of a simple function
double T[][];
T=new double[n][n];
double h[];
h=new double[n];
//vector<double> h(n,0);
for(i=0;i<n;i++)
{
  h[i]=0;
  for(int j=0;j<n;j++)
    T[i][j]=0;
}
h[0]=h0;
double r=0.5;
for( i=1;i<n;i++)
{
h[i]=h0*Math.pow(r,i);
}

for(i=0;i<n;i++)
{
T[i][0]=( P((x + h[i]),v) - P(( x - h[i]),v))/(2.0*h[i]);
}
for(m=1;m<n;m++)
{
  for(i=0;i<n-m;i++)
  {
  T[i][m]=(h[i]*h[i]*T[i+1][m-1] - h[i+m]*h[i+m]*T[i][m-1])/
          (h[i]*h[i]- h[i+m]*h[i+m]);
  }
}
double xx=T[0][n-1];
return xx;
}

public  double dPdV1(double t,double x)
{
// numerically calculated dP/dv
double h0=0.0256;
int i,m;
int n=7;
//derivative of a simple function
double T[][];
T=new double[n][n];
double h[];
h=new double[n];
//vector<double> h(n,0);
for(i=0;i<n;i++)
{
  h[i]=0;
  for(int j=0;j<n;j++)
    T[i][j]=0;
}
h[0]=h0;
double r=0.5;
for( i=1;i<n;i++)
{
h[i]=h0*Math.pow(r,i);
}

for(i=0;i<n;i++)
{
T[i][0]=( P(t,(x + h[i])) - P(t,( x - h[i])))/(2.0*h[i]);
}
for(m=1;m<n;m++)
{
  for(i=0;i<n-m;i++)
  {
  T[i][m]=(h[i]*h[i]*T[i+1][m-1] - h[i+m]*h[i+m]*T[i][m-1])/
          (h[i]*h[i]- h[i+m]*h[i+m]);
  }
}
double xx=T[0][n-1];
return xx;
}

public double dPdV(double t, double v)
{
//derivative of equation of state		
double T=t+273.15;
double Tr=T/rd.Tc;
double V=v*rd.M;
double ror=rd.Vc/V;
double pi=0;
double Tr2=Tr*Tr;double Tr3=Tr2*Tr;double Tr4=Tr2*Tr2;double Tr6=Tr2*Tr4;
double ror2=ror*ror;double ror3=ror2*ror;double ror4=ror2*ror2;double ror5=ror3*ror2;
double ror6=ror4*ror2;double ror7=ror6*ror;double ror8=ror6*ror2;double ror9=ror8*ror;
double ror10=ror8*ror2;double ror11=ror10*ror;
if(rd.refrigerantName.equals("R22"))  
{pi=
(rd.A[0]+rd.A[1]/Tr+rd.A[2]/Tr4+rd.A[3]/Tr6)+
2.0*(rd.A[4]+rd.A[5]/Tr+rd.A[6]/Tr4)*ror+
3.0*(rd.A[7]+rd.A[8]/Tr+rd.A[9]/Tr2+rd.A[10]/Tr4)*ror2+
4.0*(rd.A[11]/Tr+rd.A[12]/Tr4)*ror3+
5.0*(rd.A[13]/Tr2+rd.A[14]/Tr3+rd.A[15]/Tr4+rd.A[16]/Tr6)*ror4+
6.0*(rd.A[17])*ror5+
7.0*(rd.A[18]/Tr2+rd.A[19]/Tr6)*ror6+
8.0*(rd.A[20]/Tr)*ror7+
9.0*(rd.A[21]/Tr2+rd.A[22]/Tr6)*ror8+
10.0*(rd.A[23]/Tr+rd.A[24]/Tr6)*ror9+
11.0*(rd.A[25]/Tr+rd.A[26]/Tr2+rd.A[27]/Tr6)*ror10; 
}
else{
pi=(rd.A[0]+rd.A[1]/Tr+rd.A[2]/Tr2+rd.A[3]/Tr3)+
2.0*(rd.A[4]+rd.A[5]/Tr+rd.A[6]/Tr2)*ror+
3.0*(rd.A[7]+rd.A[8]/Tr+rd.A[9]/Tr2)*ror2+
4.0*(rd.A[10]+rd.A[11]/Tr)*ror3+
6.0*(rd.A[12]+rd.A[13]/Tr)*ror5+
8.0*(rd.A[14])*ror7;}
return -(1.0/v/v)*(P(t,v)*v+ror*pi*(rd.R/rd.M)*T);
}

// ======= v relations==========
public double v_tp(double t,double Pi)
{
//specific volume as a function of temperature and pressure	
// x in f(x)=0 given also df(x)/dx 
// Newton - Raphson type method 
double T=t+273.15;
double x=(T*rd.R/Pi)/rd.M;
int nmax=500;
double tolerance=1.0e-10;
double fx,dfx;
for(int i=0;i<nmax;i++)
{
fx=P(t,x)-Pi;
dfx=dPdV(t,x);
x-=fx/dfx;
//System.out.println("fx="+fx+"dfx="+dfx"+"x="+x);
if(Math.abs(fx)<tolerance)
{ return x; }
}
return x;
}

public double[] v_tp_braket(double t,double Pi,double x1,double x2)
{
	//returns a bracketed root 
	//if a root is not bracketed, enlarge region outwards each time
	int NTRY=50;
	double a[]=new double[2];
	double FACTOR=1.343253;
	int j;
	double f1,f2;
	if (x1 == x2) System.out.println("Bad initial range in ZBRAC");
	f1=P(t,x1)-Pi;
	f2=P(t,x2)-Pi;
	for (j=1;j<=NTRY;j++) {
		if (f1*f2 < 0.0) {a[0]=x1;a[1]=x2;return a;}
		if (Math.abs(f1) < Math.abs(f2))
			f1=P(t,(x1 += FACTOR*(x1-x2)))-Pi;
		else
			f2=P(t,(x2 += FACTOR*(x2-x1)))-Pi;
	}
	{a[0]=x1;a[1]=x2;return a;}
}

public double v_tp2(double t,double Pi)
{
//specific volume as a function of temperature and pressure	
	double T=t+273.15;
	double x=(T*rd.R/Pi)/rd.M;
	double dx=1.5436543624e-3;
	double x1=x-dx;
	double x2=x+dx;
	double aa[]=v_tp_braket(t,Pi,x1,x2);
	x1=aa[0];
	x2=aa[1];
	double ITMAX=100;
    double EPS=3.0e-8;
	int iter;
	double a=x1,b=x2,c,d,e,min1,min2;
	double tol=1.0e-5;
	c=0;
	d=0;
	e=0;
	double fa=P(t,a)-Pi,fb=P(t,b)-Pi,fc,p,q,r,s,tol1,xm;
	if (fb*fa > 0.0) System.out.println("xxx Root must be bracketed in ZBRENT");
	fc=fb;
	for (iter=1;iter<=ITMAX;iter++) {
		if (fb*fc > 0.0) {
			c=a;
			fc=fa;
			e=d=b-a;
		}
		if (Math.abs(fc) < Math.abs(fb)) {
			a=b;
			b=c;
			c=a;
			fa=fb;
			fb=fc;
			fc=fa;
		}
		tol1=2.0*EPS*Math.abs(b)+0.5*tol;
		xm=0.5*(c-b);
		if (Math.abs(xm) <= tol1 || fb == 0.0) return b;
		if (Math.abs(e) >= tol1 && Math.abs(fa) > Math.abs(fb)) {
			s=fb/fa;
			if (a == c) {
				p=2.0*xm*s;
				q=1.0-s;
			} else {
				q=fa/fc;
				r=fb/fc;
				p=s*(2.0*xm*q*(q-r)-(b-a)*(r-1.0));
				q=(q-1.0)*(r-1.0)*(s-1.0);
			}
			if (p > 0.0)  q = -q;
			p=Math.abs(p);
			min1=3.0*xm*q-Math.abs(tol1*q);
			min2=Math.abs(e*q);
			if (2.0*p < (min1 < min2 ? min1 : min2)) {
				e=d;
				d=p/q;
			} else {
				d=xm;
				e=d;
			}
		} else {
			d=xm;
			e=d;
		}
		a=b;
		fa=fb;
		if (Math.abs(d) > tol1)
			b += d;
		else
			b += (xm > 0.0 ? Math.abs(tol1) : -Math.abs(tol1));
		fb=P(t,b)-Pi;
	}
	System.out.println("Maximum number of iterations exceeded in BRENT");
return b;
}

public double[] v_th_braket(double t,double hi,double x1,double x2)
{
	//returns a bracketed root 
	//if a root is not bracketed, enlarge region outwards each time
	int NTRY=50;
	double a[]=new double[2];
	double FACTOR=1.343253;
	int j;
	double f1,f2;
	if (x1 == x2) System.out.println("Bad initial range in ZBRAC");
	f1=h(t,x1)-hi;
	f2=h(t,x2)-hi;
	for (j=1;j<=NTRY;j++) {
		if (f1*f2 < 0.0) {a[0]=x1;a[1]=x2;return a;}
		if (Math.abs(f1) < Math.abs(f2))
			f1=h(t,(x1 += FACTOR*(x1-x2)))-hi;
		else
			f2=h(t,(x2 += FACTOR*(x2-x1)))-hi;
	}
	{a[0]=x1;a[1]=x2;return a;}
}

public double v_th(double t,double hi)
{
//specific volume as a function of temperature and enthalpy	
	double T=t+273.15;
	double x=vg_app(t)+0.1;
	double dx=1.5436543624e-3;
	double x1=x-dx;
	double x2=x+dx;
	//find a bracket with the root
	double aa[]=v_th_braket(t,hi,x1,x2);
	x1=aa[0];
	x2=aa[1];
	double ITMAX=100;
    double EPS=3.0e-8;
	int iter;
	double a=x1,b=x2,c,d,e,min1,min2;
	double tol=1.0e-5;
	c=0;
	d=0;
	e=0;
	double fa=h(t,a)-hi,fb=h(t,b)-hi,fc,p,q,r,s,tol1,xm;
	if (fb*fa > 0.0) System.out.println("Root must be bracketed in ZBRENT");
	fc=fb;
	for (iter=1;iter<=ITMAX;iter++) {
		if (fb*fc > 0.0) {
			c=a;
			fc=fa;
			e=d=b-a;
		}
		if (Math.abs(fc) < Math.abs(fb)) {
			a=b;
			b=c;
			c=a;
			fa=fb;
			fb=fc;
			fc=fa;
		}
		tol1=2.0*EPS*Math.abs(b)+0.5*tol;
		xm=0.5*(c-b);
		if (Math.abs(xm) <= tol1 || fb == 0.0) return b;
		if (Math.abs(e) >= tol1 && Math.abs(fa) > Math.abs(fb)) {
			s=fb/fa;
			if (a == c) {
				p=2.0*xm*s;
				q=1.0-s;
			} else {
				q=fa/fc;
				r=fb/fc;
				p=s*(2.0*xm*q*(q-r)-(b-a)*(r-1.0));
				q=(q-1.0)*(r-1.0)*(s-1.0);
			}
			if (p > 0.0)  q = -q;
			p=Math.abs(p);
			min1=3.0*xm*q-Math.abs(tol1*q);
			min2=Math.abs(e*q);
			if (2.0*p < (min1 < min2 ? min1 : min2)) {
				e=d;
				d=p/q;
			} else {
				d=xm;
				e=d;
			}
		} else {
			d=xm;
			e=d;
		}
		a=b;
		fa=fb;
		if (Math.abs(d) > tol1)
			b += d;
		else
			b += (xm > 0.0 ? Math.abs(tol1) : -Math.abs(tol1));
		fb=h(t,b)-hi;
	}
	System.out.println("Maximum number of iterations exceeded in BRENT");
return b;
}
public double[] v_tu_braket(double t,double ui,double x1,double x2)
{
	//returns a bracketed root 
	//if a root is not bracketed, enlarge region outwards each time
	int NTRY=50;
	double a[]=new double[2];
	double FACTOR=1.343253;
	int j;
	double f1,f2;
	if (x1 == x2) System.out.println("Bad initial range in ZBRAC");
	f1=u(t,x1)-ui;
	f2=u(t,x2)-ui;
	for (j=1;j<=NTRY;j++) {
		if (f1*f2 < 0.0) {a[0]=x1;a[1]=x2;return a;}
		if (Math.abs(f1) < Math.abs(f2))
			f1=u(t,(x1 += FACTOR*(x1-x2)))-ui;
		else
			f2=u(t,(x2 += FACTOR*(x2-x1)))-ui;
	}
	{a[0]=x1;a[1]=x2;return a;}
}


public double v_tu(double t,double ui)
{
//specific volume as a function of temperature and internal energy	
	double T=t+273.15;
	double x=vg_app(t)+0.1;
	double dx=1.5436543624e-3;
	double x1=x-dx;
	double x2=x+dx;
	//find a bracket with the root
	double aa[]=v_tu_braket(t,ui,x1,x2);
	x1=aa[0];
	x2=aa[1];
	double ITMAX=100;
    double EPS=3.0e-8;
	int iter;
	double a=x1,b=x2,c,d,e,min1,min2;
	double tol=1.0e-5;
	c=0;
	d=0;
	e=0;
	double fa=u(t,a)-ui,fb=u(t,b)-ui,fc,p,q,r,s,tol1,xm;
	if (fb*fa > 0.0) System.out.println("Root must be bracketed in ZBRENT");
	fc=fb;
	for (iter=1;iter<=ITMAX;iter++) {
		if (fb*fc > 0.0) {
			c=a;
			fc=fa;
			e=d=b-a;
		}
		if (Math.abs(fc) < Math.abs(fb)) {
			a=b;
			b=c;
			c=a;
			fa=fb;
			fb=fc;
			fc=fa;
		}
		tol1=2.0*EPS*Math.abs(b)+0.5*tol;
		xm=0.5*(c-b);
		if (Math.abs(xm) <= tol1 || fb == 0.0) return b;
		if (Math.abs(e) >= tol1 && Math.abs(fa) > Math.abs(fb)) {
			s=fb/fa;
			if (a == c) {
				p=2.0*xm*s;
				q=1.0-s;
			} else {
				q=fa/fc;
				r=fb/fc;
				p=s*(2.0*xm*q*(q-r)-(b-a)*(r-1.0));
				q=(q-1.0)*(r-1.0)*(s-1.0);
			}
			if (p > 0.0)  q = -q;
			p=Math.abs(p);
			min1=3.0*xm*q-Math.abs(tol1*q);
			min2=Math.abs(e*q);
			if (2.0*p < (min1 < min2 ? min1 : min2)) {
				e=d;
				d=p/q;
			} else {
				d=xm;
				e=d;
			}
		} else {
			d=xm;
			e=d;
		}
		a=b;
		fa=fb;
		if (Math.abs(d) > tol1)
			b += d;
		else
			b += (xm > 0.0 ? Math.abs(tol1) : -Math.abs(tol1));
		fb=u(t,b)-ui;
	}
	System.out.println("Maximum number of iterations exceeded in BRENT");
return b;
}
public double[] v_ts_braket(double t,double si,double x1,double x2)
{
	//returns a bracketed root 
	//if a root is not bracketed, enlarge region outwards each time
	int NTRY=50;
	double a[]=new double[2];
	double FACTOR=1.343253;
	int j;
	double f1,f2;
	if (x1 == x2) System.out.println("Bad initial range in ZBRAC");
	f1=s(t,x1)-si;
	f2=s(t,x2)-si;
	for (j=1;j<=NTRY;j++) {
		if (f1*f2 < 0.0) {a[0]=x1;a[1]=x2;return a;}
		if (Math.abs(f1) < Math.abs(f2))
			f1=s(t,(x1 += FACTOR*(x1-x2)))-si;
		else
			f2=s(t,(x2 += FACTOR*(x2-x1)))-si;
	}
	{a[0]=x1;a[1]=x2;return a;}
}


public double v_ts(double t,double si)
{
//specific volume as a function of temperature and entropy	
	double T=t+273.15;
	double x=vg_app(t)+0.1;
	double dx=1.5436543624e-3;
	double x1=x-dx;
	double x2=x+dx;
	//find a bracket with the root
	double aa[]=v_ts_braket(t,si,x1,x2);
	x1=aa[0];
	x2=aa[1];
	double ITMAX=100;
    double EPS=3.0e-8;
	int iter;
	double a=x1,b=x2,c,d,e,min1,min2;
	double tol=1.0e-5;
	c=0;
	d=0;
	e=0;
	double fa=s(t,a)-si,fb=s(t,b)-si,fc,p,q,r,s,tol1,xm;
	if (fb*fa > 0.0) System.out.println("Root must be bracketed");
	fc=fb;
	for (iter=1;iter<=ITMAX;iter++) {
		if (fb*fc > 0.0) {
			c=a;
			fc=fa;
			e=d=b-a;
		}
		if (Math.abs(fc) < Math.abs(fb)) {
			a=b;
			b=c;
			c=a;
			fa=fb;
			fb=fc;
			fc=fa;
		}
		tol1=2.0*EPS*Math.abs(b)+0.5*tol;
		xm=0.5*(c-b);
		if (Math.abs(xm) <= tol1 || fb == 0.0) return b;
		if (Math.abs(e) >= tol1 && Math.abs(fa) > Math.abs(fb)) {
			s=fb/fa;
			if (a == c) {
				p=2.0*xm*s;
				q=1.0-s;
			} else {
				q=fa/fc;
				r=fb/fc;
				p=s*(2.0*xm*q*(q-r)-(b-a)*(r-1.0));
				q=(q-1.0)*(r-1.0)*(s-1.0);
			}
			if (p > 0.0)  q = -q;
			p=Math.abs(p);
			min1=3.0*xm*q-Math.abs(tol1*q);
			min2=Math.abs(e*q);
			if (2.0*p < (min1 < min2 ? min1 : min2)) {
				e=d;
				d=p/q;
			} else {
				d=xm;
				e=d;
			}
		} else {
			d=xm;
			e=d;
		}
		a=b;
		fa=fb;
		if (Math.abs(d) > tol1)
			b += d;
		else
			b += (xm > 0.0 ? Math.abs(tol1) : -Math.abs(tol1));
		fb=s(t,b)-si;
	}
	System.out.println("Maximum number of iterations exceeded in BRENT");
return b;
}

public double v_tx(double t,double xi)
{
//specific volume as a function of temperature and quality	
return (1-xi)*vsl_t(t)+xi*vsg_t(t);
}

public double v_ph(double pi,double hi)
{
//specific volume as a function of pressure and enthalpy	
return v_tp(t_ph(pi,hi),pi);
}

public double v_pu(double pi,double ui)
{
//specific volume as a function of pressure and enthalpy	
return v_tp(t_pu(pi,ui),pi);
}

public double v_ps(double pi,double si)
{
//specific volume as a function of pressure and enthalpy	
return v_tp(t_ps(pi,si),pi);
}
// ======= t relations==========
public double t_pv(double Pi,double v)
{
// x in f(x)=0 given also df(x)/dx 
// Newton - Raphson type method 
double V=v*rd.M;
double x=(V)*Pi/rd.R-273.15;
System.out.println("starting value ="+x);
int nmax=500;
double tolerance=1.0e-15;
double fx,dfx;
for(int i=0;i<nmax;i++)
{
fx=P(x,v)-Pi;
dfx=dPdT(x,v);
x-=fx/dfx;
if(Math.abs(fx)<tolerance)
{ return x; }
}
return x;
}

public double[] t_ph_braket(double pi,double hi,double x1,double x2)
{
	//returns a bracketed root 
	//if a root is not bracketed, enlarge region outwards each time
	int NTRY=50;
	double a[]=new double[2];
	double FACTOR=1.343253;
	int j;
	double f1,f2;
	if (x1 == x2) System.out.println("Bad initial range in ZBRAC");
	f1=h(x1,v_tp(x1,pi))-hi;
	f2=h(x2,v_tp(x2,pi))-hi;
	for (j=1;j<=NTRY;j++) {
		if (f1*f2 < 0.0) {a[0]=x1;a[1]=x2;return a;}
		if (Math.abs(f1) < Math.abs(f2))
			f1=h((x1 += FACTOR*(x1-x2)),v_tp(x1,pi))-hi;
		else
			f2=h((x2 += FACTOR*(x2-x1)),v_tp(x2,pi))-hi;
	}
	{a[0]=x1;a[1]=x2;return a;}
}


public double t_ph(double pi,double hi)
{
//specific volume as a function of pressure and enthalpy	
	double x=Ts(pi);
	double dx=2.3e-1;
	double x1=x-dx;
	double x2=x+dx;
	//find a bracket with the root
	double aa[]=t_ph_braket(pi,hi,x1,x2);
	x1=aa[0];
	x2=aa[1];
	double ITMAX=100;
    double EPS=3.0e-8;
	int iter;
	double a=x1,b=x2,c,d,e,min1,min2;
	double tol=1.0e-5;
	c=0;
	d=0;
	e=0;
	double fa=h(a,v_tp(a,pi))-hi,fb=h(b,v_tp(b,pi))-hi,fc,p,q,r,s,tol1,xm;
	if (fb*fa > 0.0) System.out.println("Root must be bracketed");
	fc=fb;
	for (iter=1;iter<=ITMAX;iter++) {
		if (fb*fc > 0.0) {
			c=a;
			fc=fa;
			e=d=b-a;
		}
		if (Math.abs(fc) < Math.abs(fb)) {
			a=b;
			b=c;
			c=a;
			fa=fb;
			fb=fc;
			fc=fa;
		}
		tol1=2.0*EPS*Math.abs(b)+0.5*tol;
		xm=0.5*(c-b);
		if (Math.abs(xm) <= tol1 || fb == 0.0) return b;
		if (Math.abs(e) >= tol1 && Math.abs(fa) > Math.abs(fb)) {
			s=fb/fa;
			if (a == c) {
				p=2.0*xm*s;
				q=1.0-s;
			} else {
				q=fa/fc;
				r=fb/fc;
				p=s*(2.0*xm*q*(q-r)-(b-a)*(r-1.0));
				q=(q-1.0)*(r-1.0)*(s-1.0);
			}
			if (p > 0.0)  q = -q;
			p=Math.abs(p);
			min1=3.0*xm*q-Math.abs(tol1*q);
			min2=Math.abs(e*q);
			if (2.0*p < (min1 < min2 ? min1 : min2)) {
				e=d;
				d=p/q;
			} else {
				d=xm;
				e=d;
			}
		} else {
			d=xm;
			e=d;
		}
		a=b;
		fa=fb;
		if (Math.abs(d) > tol1)
			b += d;
		else
			b += (xm > 0.0 ? Math.abs(tol1) : -Math.abs(tol1));
		fb=h(b,v_tp(b,pi))-hi;
	}
	System.out.println("Maximum number of iterations exceeded in BRENT");
return b;
}

public double[] t_pu_braket(double pi,double ui,double x1,double x2)
{
	//returns a bracketed root 
	//if a root is not bracketed, enlarge region outwards each time
	int NTRY=50;
	double a[]=new double[2];
	double FACTOR=1.343253;
	int j;
	double f1,f2;
	if (x1 == x2) System.out.println("Bad initial range in ZBRAC");
	f1=u(x1,v_tp(x1,pi))-ui;
	f2=u(x2,v_tp(x2,pi))-ui;
	for (j=1;j<=NTRY;j++) {
		if (f1*f2 < 0.0) {a[0]=x1;a[1]=x2;return a;}
		if (Math.abs(f1) < Math.abs(f2))
			f1=u((x1 += FACTOR*(x1-x2)),v_tp(x1,pi))-ui;
		else
			f2=u((x2 += FACTOR*(x2-x1)),v_tp(x2,pi))-ui;
	}
	{a[0]=x1;a[1]=x2;return a;}
}


public double t_pu(double pi,double ui)
{
//specific volume as a function of pressure and enthalpy	
	double x=Ts(pi);
	double dx=2.3e-1;
	double x1=x-dx;
	double x2=x+dx;
	//find a bracket with the root
	double aa[]=t_pu_braket(pi,ui,x1,x2);
	x1=aa[0];
	x2=aa[1];
	double ITMAX=100;
    double EPS=3.0e-8;
	int iter;
	double a=x1,b=x2,c,d,e,min1,min2;
	double tol=1.0e-5;
	c=0;
	d=0;
	e=0;
	double fa=u(a,v_tp(a,pi))-ui,fb=u(b,v_tp(b,pi))-ui,fc,p,q,r,s,tol1,xm;
	if (fb*fa > 0.0) System.out.println("Root must be bracketed");
	fc=fb;
	for (iter=1;iter<=ITMAX;iter++) {
		if (fb*fc > 0.0) {
			c=a;
			fc=fa;
			e=d=b-a;
		}
		if (Math.abs(fc) < Math.abs(fb)) {
			a=b;
			b=c;
			c=a;
			fa=fb;
			fb=fc;
			fc=fa;
		}
		tol1=2.0*EPS*Math.abs(b)+0.5*tol;
		xm=0.5*(c-b);
		if (Math.abs(xm) <= tol1 || fb == 0.0) return b;
		if (Math.abs(e) >= tol1 && Math.abs(fa) > Math.abs(fb)) {
			s=fb/fa;
			if (a == c) {
				p=2.0*xm*s;
				q=1.0-s;
			} else {
				q=fa/fc;
				r=fb/fc;
				p=s*(2.0*xm*q*(q-r)-(b-a)*(r-1.0));
				q=(q-1.0)*(r-1.0)*(s-1.0);
			}
			if (p > 0.0)  q = -q;
			p=Math.abs(p);
			min1=3.0*xm*q-Math.abs(tol1*q);
			min2=Math.abs(e*q);
			if (2.0*p < (min1 < min2 ? min1 : min2)) {
				e=d;
				d=p/q;
			} else {
				d=xm;
				e=d;
			}
		} else {
			d=xm;
			e=d;
		}
		a=b;
		fa=fb;
		if (Math.abs(d) > tol1)
			b += d;
		else
			b += (xm > 0.0 ? Math.abs(tol1) : -Math.abs(tol1));
		fb=u(b,v_tp(b,pi))-ui;
	}
	System.out.println("Maximum number of iterations exceeded in BRENT");
return b;
}

public double[] t_ps_braket(double pi,double si,double x1,double x2)
{
	//returns a bracketed root 
	//if a root is not bracketed, enlarge region outwards each time
	int NTRY=50;
	double a[]=new double[2];
	double FACTOR=1.343253;
	int j;
	double f1,f2;
	if (x1 == x2) System.out.println("Bad initial range in ZBRAC");
	f1=s(x1,v_tp(x1,pi))-si;
	f2=s(x2,v_tp(x2,pi))-si;
	for (j=1;j<=NTRY;j++) {
		if (f1*f2 < 0.0) {a[0]=x1;a[1]=x2;return a;}
		if (Math.abs(f1) < Math.abs(f2))
			f1=s((x1 += FACTOR*(x1-x2)),v_tp(x1,pi))-si;
		else
			f2=s((x2 += FACTOR*(x2-x1)),v_tp(x2,pi))-si;
	}
	{a[0]=x1;a[1]=x2;return a;}
}

public double t_ps(double pi,double si)
{
//specific volume as a function of pressure and enthalpy	
	double x=Ts(pi);
	double dx=2.3e-1;
	double x1=x-dx;
	double x2=x+dx;
	//find a bracket with the root
	double aa[]=t_ps_braket(pi,si,x1,x2);
	x1=aa[0];
	x2=aa[1];
	double ITMAX=100;
    double EPS=3.0e-8;
	int iter;
	double a=x1,b=x2,c,d,e,min1,min2;
	double tol=1.0e-5;
	c=0;
	d=0;
	e=0;
	double fa=s(a,v_tp(a,pi))-si,fb=s(b,v_tp(b,pi))-si,fc,p,q,r,s,tol1,xm;
	if (fb*fa > 0.0) System.out.println("Root must be bracketed");
	fc=fb;
	for (iter=1;iter<=ITMAX;iter++) {
		if (fb*fc > 0.0) {
			c=a;
			fc=fa;
			e=d=b-a;
		}
		if (Math.abs(fc) < Math.abs(fb)) {
			a=b;
			b=c;
			c=a;
			fa=fb;
			fb=fc;
			fc=fa;
		}
		tol1=2.0*EPS*Math.abs(b)+0.5*tol;
		xm=0.5*(c-b);
		if (Math.abs(xm) <= tol1 || fb == 0.0) return b;
		if (Math.abs(e) >= tol1 && Math.abs(fa) > Math.abs(fb)) {
			s=fb/fa;
			if (a == c) {
				p=2.0*xm*s;
				q=1.0-s;
			} else {
				q=fa/fc;
				r=fb/fc;
				p=s*(2.0*xm*q*(q-r)-(b-a)*(r-1.0));
				q=(q-1.0)*(r-1.0)*(s-1.0);
			}
			if (p > 0.0)  q = -q;
			p=Math.abs(p);
			min1=3.0*xm*q-Math.abs(tol1*q);
			min2=Math.abs(e*q);
			if (2.0*p < (min1 < min2 ? min1 : min2)) {
				e=d;
				d=p/q;
			} else {
				d=xm;
				e=d;
			}
		} else {
			d=xm;
			e=d;
		}
		a=b;
		fa=fb;
		if (Math.abs(d) > tol1)
			b += d;
		else
			b += (xm > 0.0 ? Math.abs(tol1) : -Math.abs(tol1));
		fb=s(b,v_tp(b,pi))-si;
	}
	System.out.println("Maximum number of iterations exceeded in BRENT");
return b;
}


// ======= h relations =============

public double h1(double t, double v)
{
double T=t+273.15;
double Tr=T/rd.Tc;
double Tc=rd.Tc;double Tc2=Tc*Tc;double Tc3=Tc2*Tc;double Tc4=Tc2*Tc2;
double V=v*rd.M;
double ror=rd.Vc/V;
double R=rd.R/rd.M;
double pi=0;
double pi1=0;
double pi2=0;
double T0=273.15/Tc;
double T02=T0*T0;double T03=T0*T02;double T04=T02*T02;double T05=T04*T0;
double Tr2=Tr*Tr;double Tr3=Tr2*Tr;double Tr4=Tr2*Tr2;double Tr5=Tr*Tr4;double Tr6=Tr2*Tr4;
double ror2=ror*ror;double ror3=ror2*ror;double ror4=ror2*ror2;double ror5=ror3*ror2;
double ror6=ror4*ror2;double ror7=ror6*ror;double ror8=ror6*ror2;double ror9=ror8*ror;
double ror10=ror8*ror2;double ror11=ror10*ror;
if(rd.refrigerantName.equals("R22"))  
{pi=(rd.A[0]+2.0*rd.A[1]/Tr+5.0*rd.A[2]/Tr4+7.0*rd.A[3]/Tr6)*ror;
 pi+=(2.0*rd.A[4]+3.0*rd.A[5]/Tr+6.0*rd.A[6]/Tr4)*ror2/2.0;
 pi+=(3.0*rd.A[7]+4.0*rd.A[8]/Tr+5.0*rd.A[9]/Tr2+7.0*rd.A[10]/Tr4)*ror3/3.0;
 pi+=(5.0*rd.A[11]/Tr+8.0*rd.A[12]/Tr4)*ror4/4.0;
 pi+=(7.0*rd.A[13]/Tr2+8.0*rd.A[14]/Tr3+9.0*rd.A[15]/Tr4+11.0*rd.A[16]/Tr6)*ror5/5.0;
 pi+=(rd.A[17])*ror6;
 pi+=(9.0*rd.A[18]/Tr2+13.0*rd.A[19]/Tr6)*ror7/7.0;
 pi+=(9.0*rd.A[20]/Tr)*ror8/8.0;
 pi+=(11.0*rd.A[21]/Tr2+15.0*rd.A[22]/Tr6)*ror9/9.0;
 pi+=(11.0*rd.A[23]/Tr+16.0*rd.A[24]/Tr6)*ror10/10.0;
 pi+=(12.0*rd.A[25]/Tr+13.0*rd.A[26]/Tr2+17.0*rd.A[27]/Tr6)*ror11/11.0;
 pi*=R*T;
}
else{
pi=R*T*((rd.A[0]+2.0*rd.A[1]/Tr+3.0*rd.A[2]/Tr2+4.0*rd.A[3]/Tr3)*ror+
(2.0*rd.A[4]+3.0*rd.A[5]/Tr+4.0*rd.A[6]/Tr2)*ror2/2.0+
(3.0*rd.A[7]+4.0*rd.A[8]/Tr+5.0*rd.A[9]/Tr2)*ror3/3.0+
(4.0*rd.A[10]+5.0*rd.A[11]/Tr)*ror4/4.0+
(6.0*rd.A[12]+7.0*rd.A[13]/Tr)*ror6/6.0+
(rd.A[14])*ror8)+dH0(t)+rd.X;}
pi=pi+dH0(t)+rd.X;
return pi;
}

public double h(double t, double v)
{
double h1=h1(t,v);
double x=x_tv(t,v);
if(x>1)
{
 h1=h1(t,v);
}
else if(x<0)
{
 double P=P(t,v);
 h1=hsl_t(t)+vsl_t(t)*(P-Ps(t));
 //h1=h1(t,v)-hslg_t(t);
}
else if(x==0)
{
 h1=hsl_t(t);
}
else if(x==1)
{
 h1=hsg_t(t);
}
else if(x>0 || x<1)
{	
 h1=hsl_t(t)*(1-x)+hsg_t(t)*x;
}
return h1;
}

// ======= u relations =============
public double u(double t, double v)
{return h(t,v)-P(t,v)*v;}

// ======= s relations =============
public double s1(double t, double v)
{
// t degree C	
// v m^3/kg
// s1 kJ/kg K
double T=t+273.15;
double Tr=T/rd.Tc;
double Tc=rd.Tc;
double V=v*rd.M;
double ror=rd.Vc/V;
double ds=0;
double Tr2=Tr*Tr;double Tr3=Tr2*Tr;double Tr4=Tr2*Tr2;double Tr5=Tr4*Tr;
double Tr6=Tr2*Tr4;double Tr7=Tr6*Tr;
double Tc2=Tc*Tc;double Tc3=Tc2*Tc;double Tc4=Tc2*Tc2;
double ror2=ror*ror;double ror3=ror2*ror;double ror4=ror2*ror2;double ror5=ror3*ror2;
double ror6=ror4*ror2;double ror7=ror6*ror;double ror8=ror6*ror2;double ror9=ror8*ror;
double ror10=ror8*ror2;double ror11=ror10*ror;
double T0=273.15/Tc;
double T02=T0*T0;double T03=T0*T02;double T04=T02*T02;double T05=T04*T0;
double B[]=new double[11];
double dB[]=new double[11];
double R=rd.R/rd.M;
if(rd.refrigerantName.equals("R22"))  
{B[0]=(rd.A[0]+rd.A[1]/Tr+rd.A[2]/Tr4+rd.A[3]/Tr6);
 B[1]=(rd.A[4]+rd.A[5]/Tr+rd.A[6]/Tr4);
 B[2]=(rd.A[7]+rd.A[8]/Tr+rd.A[9]/Tr2+rd.A[10]/Tr4);
 B[3]=(rd.A[11]/Tr+rd.A[12]/Tr4);
 B[4]=(rd.A[13]/Tr2+rd.A[14]/Tr3+rd.A[15]/Tr4+rd.A[16]/Tr6);
 B[5]=(rd.A[17]);
 B[6]=(rd.A[18]/Tr2+rd.A[19]/Tr6);
 B[7]=(rd.A[20]/Tr);
 B[8]=(rd.A[21]/Tr2+rd.A[22]/Tr6);
 B[9]=(rd.A[23]/Tr+rd.A[24]/Tr6);
 B[10]=(rd.A[25]/Tr+rd.A[26]/Tr2+rd.A[27]/Tr6);
 dB[0]=(-rd.A[1]/Tr2-4.0*rd.A[2]/Tr5-6.0*rd.A[3]/Tr7);
 dB[1]=(-rd.A[5]/Tr2-4.0*rd.A[6]/Tr5);
 dB[2]=(-rd.A[8]/Tr2-2.0*rd.A[9]/Tr3-4.0*rd.A[10]/Tr5);
 dB[3]=(-rd.A[11]/Tr2-4.0*rd.A[12]/Tr5);
 dB[4]=(-2.0*rd.A[13]/Tr3-3.0*rd.A[14]/Tr4-4.0*rd.A[15]/Tr5-6.0*rd.A[16]/Tr7);
 dB[5]=0.0;
 dB[6]=(-2.0*rd.A[18]/Tr3-6.0*rd.A[19]/Tr7);
 dB[7]=(-rd.A[20]/Tr2);
 dB[8]=(-2.0*rd.A[21]/Tr3-6.0*rd.A[22]/Tr7);
 dB[9]=(-rd.A[23]/Tr2-6.0*rd.A[24]/Tr7);
 dB[10]=(-rd.A[25]/Tr2-2.0*rd.A[26]/Tr3-6.0*rd.A[27]/Tr7);	
}
else{
 B[0]=(rd.A[0]+rd.A[1]/Tr+rd.A[2]/Tr2+rd.A[3]/Tr3);
 B[1]=(rd.A[4]+rd.A[5]/Tr+rd.A[6]/Tr2);
 B[2]=(rd.A[7]+rd.A[8]/Tr+rd.A[9]/Tr2);
 B[3]=(rd.A[10]+rd.A[11]/Tr);
 B[4]=0.0;
 B[5]=(rd.A[12]+rd.A[13]/Tr);
 B[6]=0.0;
 B[7]=(rd.A[14]);
 B[8]=0.0;
 B[9]=0.0;
 B[10]=0.0;	
 dB[0]=(-rd.A[1]/Tr2-2.0*rd.A[2]/Tr3-3.0*rd.A[3]/Tr4);
 dB[1]=(-rd.A[5]/Tr2-2.0*rd.A[6]/Tr3);
 dB[2]=(-rd.A[8]/Tr2-2.0*rd.A[9]/Tr3);
 dB[3]=(-rd.A[11]/Tr2);
 dB[4]=0.0;
 dB[5]=(-rd.A[13]/Tr2);
 dB[6]=0.0;
 dB[7]=0.0;
 dB[8]=0.0;
 dB[9]=0.0;
 dB[10]=0.0;
}
ds=Math.log(ror);
for(int i=0;i<11;i++)
{
ds+=(B[i]+Tr*dB[i])*Math.pow(ror,(i+1))/(i+1);
}
ds*=-R;
ds+=dS0(t)+rd.Y;

/*
if(rd.refrigerantName.equals("R22"))  
{pi=R*(-Math.log(ror)+
(rd.A[0]-3.0*rd.A[2]/Tr4-5.0*rd.A[3]/Tr6)*ror+
(rd.A[4]-3.0*rd.A[6]/Tr2)*ror2/2.0+
(rd.A[7]-rd.A[9]/Tr2-3.0*rd.A[10]/Tr4)*ror3/3.0+
(-3.0*rd.A[12]/Tr4)*ror4/4.0 -
(rd.A[13]/Tr2+2.0*rd.A[14]/Tr3+3.0*rd.A[15]/Tr4+5.0*rd.A[16]/Tr6)*ror5/5+
(rd.A[17])*ror6/6.0-
(rd.A[18]/Tr2+5.0*rd.A[19]/Tr6)*ror7/7.0-
(rd.A[21]/Tr2+5.0*rd.A[22]/Tr6)*ror9/9.0-
(rd.A[24]/Tr6)*ror10/2.0-
(rd.A[26]/Tr2+5.0*rd.A[27]/Tr6)*ror11/11.0);
pi=pi+dS0(t)+rd.Y;
//Math.log(Tr/T0))-
//((rd.G[1]-R)*Math.log(Tr/T0)/2.0+rd.G[2]*Tc*(Tr-T0)+rd.G[3]*Tc2*(Tr2-T02)/2.0+
//rd.G[4]*Tc3*(Tr3-T03)/3.0+rd.G[5]*Tc4*(Tr4-T04)/4.0)+rd.Y; 
}
else{
pi=-R*(Math.log(ror)+(rd.A[0]-rd.A[2]/Tr2-2.0*rd.A[3]/Tr3)*ror+
(rd.A[4]-rd.A[6]/Tr2)*ror2/2.0+
(rd.A[7]-rd.A[9]/Tr2)*ror3/3.0+
(rd.A[10])*ror4/4.0+
(rd.A[12])*ror6/6.0)+
+R*((rd.A[15])*ror8+Tc*((rd.G[0]/Tc)*Math.log(Tr/T0)+
rd.G[2]*Tc*(Tr-T0)+rd.G[3]*Tc3*(Tr3-T03)/3.0+
rd.G[5]*Tc4*(Tr4-T04)/4.0))+rd.Y;}
return pi;
*/
return ds;
}

public double s(double t, double v)
{
double s1=s1(t,v);
double x=x_tv(t,v);
if(x>1)
{
 s1=s1(t,v);
}
else if(x<0)
{
 s1=s1(t,v)-sslg_t(t);
}
else if(x==0)
{
 s1=ssl_t(t);
}
else if(x==1)
{
 s1=ssg_t(t);
}
else if(x>0 || x<1)
{	
 s1=ssl_t(t)*(1-x)+ssg_t(t)*x;
}
return s1;
}

// ======= ideal gas relations =============
public double Cv(double t)
{//ideal gas heat capacity
return Cp(t)-rd.R/rd.M; //KJ/kmol K
}

public double Cp(double t)
{//ideal gas heat capacity
double T=t+273.15;
return rd.G[0]/T+rd.G[1]+rd.G[2]*T+rd.G[3]*T*T+rd.G[4]*T*T*T+rd.G[5]*T*T*T*T; 
}

public double dH0(double t)
{//ideal gas enthalpy difference
double T=t+273.15;
double T0=273.15;
return rd.G[0]*Math.log(T/T0)+rd.G[1]*(T-T0)+rd.G[2]/2.0*(T*T-T0*T0)+
rd.G[3]/3.0*(T*T*T-T0*T0*T0)+rd.G[4]/4.0*(T*T*T*T-T0*T0*T0*T0)+
rd.G[5]/5.0*(T*T*T*T*T-T0*T0*T0*T0*T0); 
}

public double dS0(double t)
{//ideal gas entrophy difference
double T=t+273.15;
double T0=273.15;
double R=rd.R/rd.M;
return -rd.G[0]*(1.0/T0-1.0/T)+(rd.G[1]-R)*Math.log(T/T0)+rd.G[2]*(T-T0)+rd.G[3]/2.0*(T*T-T0*T0)+
rd.G[4]/3.0*(T*T*T-T0*T0*T0)+rd.G[5]/4.0*(T*T*T*T-T0*T0*T0*T0); 
}
// saturation state
//=========================================================
public String phase(double x)
{
//phase of the steam-water system
if( x==3)
  return "kritik üstü sıvı";
else
 {
 if(x<0)
   return "sıvı";
 else if( x==0.0)
   return "doymuş sıvı";
 else if(x==1.0)
   return "doymuş buhar";
 else if( x>0.0 && x<1.0)
   return "doymuş sıvı-buhar karışımı";
 else // (x>1)
   return "kızgın buhar";
 }
}
public String phase(double t, double v)
{
//phase of the vapour-liquid system
if( Ps(t) > rd.Pc)
  return "kritik üstü sıvı";
else
 {
 double ts=t;
 double vsl=vsl_t(ts);
 double vsg=vsg_t(ts);
 if(v<vsl)
   return "sıvı";
 else if(v==vsl)
   return "doymuş sıvı";
 else if(v==vsg)
   return "doymuş buhar";
 else if(v>vsl && v<vsg)
   return "doymuş sıvı-buhar karışımı";
 else // (v>vsg)
   return "kızgın buhar";
 }
}

public double x_tv(double t, double v)
{
//phase and/or quality of the steam-water system
double xx=0;
double p=P(t,v);
double Ps=Ps(t);
double vsl=vsl_t(t);
double vsg=vsg_t(t);
if(Ps >= rd.Pc*1e3)
  xx=3;
if(p < Ps)  xx=2;
else if (p > Ps) xx=-1;
else //p==Ps
{
  if(v==vsl)
    xx=0;
  else if(v==vsg)
    xx=1;
  else if(v>vsl && v<vsg)
    xx=(v-vsl)/(vsg-vsl);
}    
  return xx;  
}

public double x_pv(double p, double v)
{
//phase and/or quality of the two pahase system
double xx=0;
double t=t_pv(p,v);
double Ps=Ps(t);
double vsl=vsl_t(t);
double vsg=vsg_t(t);
if(Ps >= rd.Pc*1e3)
  xx=3;
if(p < Ps)  xx=2;
else if (p > Ps) xx=-1;
else //p==Ps
{
  if(v==vsl)
    xx=0;
  else if(v==vsg)
    xx=1;
  else if(v>vsl && v<vsg)
    xx=(v-vsl)/(vsg-vsl);
}    
return xx;  
}

      public double Ps(double t)
      {
      // saturation pressure as a function of saturation temperature
      double T=t+273.15;
      if( T >= rd.Tc ) return rd.Pc;
      double Tr=T/rd.Tc;
      double to=1-Tr;
      double to3=to*to*to;double to4=to3*to;double to5=to4*to;
      double pi=1.0/Tr*(rd.aVP[0]*to+rd.aVP[1]*Math.pow(to,1.5)+rd.aVP[2]*to3+rd.aVP[3]*to4+rd.aVP[4]*to5);
      pi=rd.Pc*1e3*Math.exp(pi);
      return pi;
      }

      public double dPsdT(double t)
      {
      //derivative of saturation pressure as a function of saturation temperature
      double T=t+273.15;
      double Tr=T/rd.Tc;
      double to=1-Tr;
      double tom1=to-1.0;double tom1_2=tom1*tom1;
      double to2=to*to;double to3=to*to*to;double to4=to3*to;double to5=to4*to;
      double pi=1.0/rd.Tc*(1/tom1*(rd.aVP[0]+1.5*rd.aVP[1]*Math.sqrt(to)+3.0*rd.aVP[2]*to2+4.0*rd.aVP[3]*to3+5.0*rd.aVP[4]*to4)+
      -1.0/tom1_2*(rd.aVP[0]*to+rd.aVP[1]*Math.pow(to,1.5)+rd.aVP[2]*to3+rd.aVP[3]*to4+rd.aVP[4]*to5));
      pi=Ps(t)*pi;
      return pi; //MPa
      }
      
public double Ts(double Pi)
{
// x in f(x)=0 given also df(x)/dx 
// Newton - Raphson type method 
if( Pi >= rd.Pc ) return rd.Tc;
double x=Ts1(Pi);
//double x=(rd.Tc-273.15);
int nmax=500;
double tolerance=1.0e-15;
double fx,dfx;
for(int i=0;i<nmax;i++)
{
fx=Ps(x)-Pi;
dfx=dPsdT(x);
//System.out.println("x="+x+"fx="+fx+"dfx="+dfx+"dfx1="+dfx1);
x-=fx/dfx;
if(Math.abs(fx)<tolerance)
{ return x; }
}
return x;
}

public  double dPsdT1(double x)
{
// numeric dP/dv
double h0=0.0256;
int i,m;
int n=7;
//derivative of a simple function
double T[][];
T=new double[n][n];
double h[];
h=new double[n];
//vector<double> h(n,0);
for(i=0;i<n;i++)
{
  h[i]=0;
  for(int j=0;j<n;j++)
    T[i][j]=0;
}
h[0]=h0;
double r=0.5;
for( i=1;i<n;i++)
{
h[i]=h0*Math.pow(r,i);
}

for(i=0;i<n;i++)
{
T[i][0]=( Ps(x + h[i]) - Ps(x - h[i]))/(2.0*h[i]);
}
for(m=1;m<n;m++)
{
  for(i=0;i<n-m;i++)
  {
  T[i][m]=(h[i]*h[i]*T[i+1][m-1] - h[i+m]*h[i+m]*T[i][m-1])/
          (h[i]*h[i]- h[i+m]*h[i+m]);
  }
}
double xx=T[0][n-1];
return xx;
}

public double Ts1(double Pi)
{
    // utilises brent algorithms and bracketing 
    // techniques to find the root
    double x1,x2,xl,xg;
	x1=-200;
	x2=rd.Tc-273.15+200;
	double aa[]=new double[2];
    double ITMAX=100;
    double EPS=3.0e-3;
	int iter;
	double a=x1,b=x2,c,d,e,min1,min2;
	double tol=1.0e-5;
	c=0;
	d=0;
	e=0;
	double fa=Ps(a)-Pi;
	double fb=Ps(b)-Pi;
	double fc,p,q,r,s,tol1,xm;
	if (fb*fa > 0.0)  System.out.println("Root must be bracketed in ZBRENT");
	fc=fb;
	for (iter=1;iter<=ITMAX;iter++) {
		if (fb*fc > 0.0) {
			c=a;
			fc=fa;
			e=d=b-a;
		}
		if (Math.abs(fc) < Math.abs(fb)) {
			a=b;
			b=c;
			c=a;
			fa=fb;
			fb=fc;
			fc=fa;
		}
		tol1=2.0*EPS*Math.abs(b)+0.5*tol;
		xm=0.5*(c-b);
		if (Math.abs(xm) <= tol1 || fb == 0.0) return b;
		if (Math.abs(e) >= tol1 && Math.abs(fa) > Math.abs(fb)) {
			s=fb/fa;
			if (a == c) {
				p=2.0*xm*s;
				q=1.0-s;
			} else {
				q=fa/fc;
				r=fb/fc;
				p=s*(2.0*xm*q*(q-r)-(b-a)*(r-1.0));
				q=(q-1.0)*(r-1.0)*(s-1.0);
			}
			if (p > 0.0)  q = -q;
			p=Math.abs(p);
			min1=3.0*xm*q-Math.abs(tol1*q);
			min2=Math.abs(e*q);
			if (2.0*p < (min1 < min2 ? min1 : min2)) {
				e=d;
				d=p/q;
			} else {
				d=xm;
				e=d;
			}
		} else {
			d=xm;
			e=d;
		}
		a=b;
		fa=fb;
		if (Math.abs(d) > tol1)
			b += d;
		else
			b += (xm > 0.0 ? Math.abs(tol1) : -Math.abs(tol1));
		fb=Ps(b)-Pi;
	}
	System.out.println("Maximum number of iterations exceeded in algorithm, root might not be valid");
return b;
}

public double rol(double t)
{//saturated liquid density
double T=t+273.15;
if( T >= rd.Tc) return rd.M/rd.Vc;
double roli=0;
roli=1.0;
double to=(rd.Tc-T)/rd.Tc;
for(int i=0;i<6;i++)
{roli+=rd.E[i]*Math.pow(to,((double)(i+1)/3.0));}
return roli/(rd.Vc/rd.M);
}

public double vis_liq(double t)
{
//viscosity of liquid Pa.s = kg/(m.s)
double visg=rd.vl[8];
  for(int i=7;i>=0;i--)
   { visg=visg*t+rd.vl[i]; }
return visg;
}

public double vis_vap(double t)
{
//viscosity of vapor Pa.s = kg/(m.s)	
double visg=rd.vv[8];
  for(int i=7;i>=0;i--)
   { visg=visg*t+rd.vv[i]; }
return visg;
}
public double k_liq(double t)
{
//thermal conductivity of liquid W/(m K)	
double kg=rd.kl[8];
  for(int i=7;i>=0;i--)
   { kg=kg*t+rd.kl[i]; }
return kg;
}

public double k_vap(double t)
{
//thermal conductivity of vapor	W/(m K)
double kg=rd.kv[8];
  for(int i=7;i>=0;i--)
   { kg=kg*t+rd.kv[i]; }
return kg;
}

public double st(double t)
{
//surface tension Pa.m = N/m	
double stg=rd.st[8];
  for(int i=7;i>=0;i--)
   { stg=stg*t+rd.st[i]; }
return stg;
}

public double alpha_liq(double t)
{
//alpha=k/(ro*Cp)
return k_liq(t)*vsl_t(t)/(Cpsl_t(t)*1e3);
}

public double alpha_vap(double t)
{
//alpha=k/(ro*Cp)
return k_vap(t)*vsg_t(t)/(Cpsg_t(t)*1e3);
}

public double Pr_liq(double t)
{
//Doymuş sıvı prandtl sayısı Pr=vis/(k*Cp)
return Cpsl_t(t)*1e3*vis_liq(t)/k_liq(t);
}
public double Pr_vap(double t)
{
//Doymuş buhar prandtl sayısı Pr=vis/(k*Cp)
return Cpsg_t(t)*1e3*vis_vap(t)/k_vap(t);
}

public double rog_app(double t)
{
if( (t+273.15) >= rd.Tc) return rd.M/rd.Vc;
//saturated vapor density (approximation)
double T=t+273.15;
double roc=1.0/rd.Vc;
double to=1-T/rd.Ti;
double ro=0;
if(   rd.refrigerantName.equals("R22"))
{  for(int i=0;i<10;i++)
       {  ro+=rd.F[i]*Math.pow(to,(double)(i)/3.0); }
   for(int i=10;i<19;i++)
       {  ro+=rd.F[i]*Math.pow(to,(double)(9-i)/3.0); }
   ro+=rd.F[19]*Math.log(T/rd.Ti);
   ro=Math.exp(ro)*rd.roi*rd.M;}
else if(rd.refrigerantName.equals("R12"))
{
   for(int i=0;i<24;i++)
       {  ro+=rd.F[i]*Math.pow(to,(double)(i-10)/3.0); }
   ro=Math.exp(ro)*rd.roi*rd.M;
}
return ro;
}

public double vg_app(double t)
{
//saturated vapor specific volume (approximation)	
return 1.0/rog_app(t);
}
public double rog(double t) {return 1.0/vsg_t(t);}

public double vsg_t(double t)
{
if( (t+273.15) >= rd.Tc) return rd.Vc/rd.M;
// x in f(x)=0 given also df(x)/dx 
// Newton - Raphson type method 
//return v_tp(t,Ps(t));
double Pi=Ps(t);
double T=t+273.15;
double x=x=vg_app(t);
int nmax=100;
double tolerance=1.0e-2;
double fx,dfx;
for(int i=0;i<nmax;i++)
{
double pp=P(t,x);	
fx=P(t,x)-Pi;
dfx=dPdV(t,x);
//System.out.println("pp="+pp+"pi="+Pi+"fx="+fx+"dfx="+dfx);
x-=fx/dfx;
if(Math.abs(fx)<tolerance)
{ return x; }
}
return x;

}

public double vsg_p(double p)
{return vsg_t(Ts(p));}

public double vg1(double t)
{
if( (t+273.15) >= rd.Tc) return rd.M/rd.Vc;
double Pi=Ps(t);
double x=vg_app(t);
int nmax=100;
double dx=1.0e-2;
double dx1=1.9565e-2;
double x1m;
double tolerance=1.0e-8;
double fx,fxm,fxm1,dfx,dfxm,d2fx;
for(int i=0;i<nmax;i++)
{
fx=P(t,x)-Pi;
fxm=P(t,(x-dx))-Pi;
fxm1=P(t,(x-dx1))-Pi;
dfx=(fx-fxm)/dx;
dfxm=(fx-fxm1)/dx1;
d2fx=-6.0/dx/dx*(fx-fxm)+2.0/dx*(2.0*dfx+dfxm);
x-=(fx/dfx+.5*fx*fx/(dfx*dfx*dfx)*d2fx);
if(Math.abs(fx)<tolerance)
{
return x;
}
}
return x;
}

public double vsl_t(double t)
{ 	return 1.0/rol(t);}

public double hslg_t1(double t)
{   if( (t+273.15) >= rd.Tc) return 0.0;
	return (t+273.13)*dPsdT1(t)*(vsg_t(t)-vsl_t(t));}
	
public double hslg_t(double t)
{   if( (t+273.15) >= rd.Tc) return 0.0;
	return (t+273.13)*dPsdT(t)*(vsg_t(t)-vsl_t(t));}

public double hsl_t(double t)
{   if( (t+273.15) >= rd.Tc) return hsg_t(rd.Tc-273.15);
	return hsg_t(t)-hslg_t(t);}

public double Cpsl_t(double t)
{ // Specific heat of saturated liquid
	return Math.abs(hsl_t(t)-hsl_t(t-0.1))/0.1;
}	
public double Cpsg_t(double t)
{ // Specific heat of saturated liquid
	return Math.abs(hsg_t(t)-hsg_t(t-0.1))/0.1;
}
public double hsg_t(double t)
{if( (t+273.15) >= rd.Tc){t=rd.Tc-273.15;};return h1(t,vsg_t(t));}

public double ssg_t(double t)
{if( (t+273.15) >= rd.Tc){t=rd.Tc-273.15;}
 return s1(t,vsg_t(t));}

public double ssl_t(double t)
{ if( (t+273.15) >= rd.Tc) return ssg_t(rd.Tc-273.15);
  return ssg_t(t)-sslg_t(t);}

public double sslg_t(double t)
{if( (t+273.15) >= rd.Tc) return 0.0;
 return hslg_t(t)/(t+273.15);}
//==========toString output methods
public String[][] toString(String s, double v1, double v2)
{
  String s1[][]=new String[9][3];
   s1[0][0]="P, basınç        ";
   s1[1][0]="T, sıcaklık      ";
   s1[2][0]="v, özgül hacim   ";
   s1[3][0]="h, entalpi       ";
   s1[4][0]="u, iç enerji     ";
   s1[5][0]="s, entropi       ";
   s1[6][0]="x, kuruluk derec.";
   s1[7][0]="   yoğunluk";
   s1[8][0]="faz              ";
   if(unit.equals("SI"))
   {
   s1[0][2]=" kPa           ";
   s1[1][2]=" derece C      ";
   s1[2][2]=" m^3/kg        ";
   s1[3][2]=" KJ/kg         ";
   s1[4][2]=" KJ/kg         ";
   s1[5][2]=" KJ/kg         ";
   s1[6][2]=" kg buh/kg kar ";
   s1[7][2]=" kg/m^3        ";
   s1[8][2]=" ";
   }
   else
   {
   s1[0][2]=" lbf/in^2       ";
   s1[1][2]=" derece F          ";
   s1[2][2]=" ft^3/lbm       ";
   s1[3][2]=" BTU/lbm        ";
   s1[4][2]=" BTU/lbm        ";
   s1[5][2]=" BTU/lbm        ";
   s1[6][2]=" lbm vap/lbm mix";
   s1[7][2]=" lbm/ft^3        ";
   s1[8][2]=" ";
   }
double pp[]=property(s,v1,v2);
for(int i=0;i<7;i++)
   {s1[i][1]=""+pp[i];}
    s1[7][1]=""+pp[7];
    s1[8][1]=phase(pp[6]);
return s1;
}

public String toString1(String s, double v1, double v2)
{
   String[] s1=new String[7];
   String[] s2=new String[7];
   String s3="";
   double pp[];
   s1[0]="P :";
   s1[1]="T :";
   s1[2]="v :";
   s1[3]="h :";
   s1[4]="u :";
   s1[5]="s :";
   s1[6]="x :";
   s1[7]="  :";
   if(unit.equals("SI"))
   {
   s2[0]=" kPa           ";
   s2[1]=" derece C      ";
   s2[2]=" m^3/kg        ";
   s2[3]=" KJ/kg         ";
   s2[4]=" KJ/kg         ";
   s2[5]=" KJ/kg         ";
   s2[6]=" kg buh/kg kar ";
   s2[7]=" kg/m^3        ";
   }
   else
   {
   s2[0]=" lbf/in^2       ";
   s2[1]=" derece F       ";
   s2[2]=" ft^3/lbm       ";
   s2[3]=" BTU/lbm        ";
   s2[4]=" BTU/lbm        ";
   s2[5]=" BTU/lbm        ";
   s2[6]=" lbm buh/lbm kar";
   s2[7]=" lbm/ft^3        ";
   }
pp=property(s,v1,v2);
s3+=phase(pp[6])+"\n";

for(int i=1;i<9;i++)
{
s3+=s1[i-1]+s2[i-1]+pp[i-1]+"\n";
}
return s3;
}

public String toString()
{
   String[] s1=new String[7];
   String[] s2=new String[7];
   String s3="";
   double pp[];
   s1[0]="P :";
   s1[1]="T :";
   s1[2]="v :";
   s1[3]="h :";
   s1[4]="u :";
   s1[5]="s :";
   s1[6]="x :";
   s1[7]="  :"; 
   if(unit.equals("SI"))
   {
   s2[0]=" kPa           ";
   s2[1]=" derece C      ";
   s2[2]=" m^3/kg        ";
   s2[3]=" KJ/kg         ";
   s2[4]=" KJ/kg         ";
   s2[5]=" KJ/kg K       ";
   s2[6]=" kg buh/kg kar ";
   s2[7]=" kg/m^3        ";
   }
   else
   {
   s2[0]=" lbf/in^2       ";
   s2[1]=" derece F       ";
   s2[2]=" ft^3/lbm       ";
   s2[3]=" BTU/lbm        ";
   s2[4]=" BTU/lbm        ";
   s2[5]=" BTU/lbm        ";
   s2[6]=" lbm buh/lbm kar";
   s2[7]=" lbm/ft^3       ";
   }
for(int i=0;i<8;i++)
{
s3+=s1[i]+s2[i]+"\n";
}
return s3;
}

void setUnit(String x)
{
//sets SI or English units
if(x.equals("EN")) unit=u[1];
else unit=u[0];
}

void setUnit(boolean x)
{
//sets SI or English units
if(x) unit=u[0];else unit=u[1];
}


//========== reading refrigerant data =============
public String readsogut(String rName)
{
String s="";
Vector gi=new Vector();

rd=new sogut1Data("R22","chlorodifluoromethane",
8.31434,5.9328,369.3,
4.988,0.168556,369.3,86.469,414.921163634314,1.459802785);
rd.EOS_in(0.54576200,  -1.3919800   ,-0.4325620,0.02214000,-0.1307418,
          0.79211000,  -0.1670240   ,0.56743874,-1.3507100,-0.1154870,
          1.02456700,   0.34435035  ,-0.4082677,0.08300990,-0.1899033,
          0.08821727,   0.01905950  ,-0.0376347,0.03329212,-0.03794234,
          0.7869090e-2,-0.4626965e-2,0.2336405e-1,-0.2066556e-2,-0.1050183e-1,
          0.5276995e-3,0.2095470e-3,0.1346363e-2);

rd.vapor_pressure_in(-7.0340913,1.4030736,-4.9605880,8.8828089,-10.6006380);
rd.ideal_gas_heat_capacity_in(0.0,0.337055,0.283027e-3,0.537315e-5,-0.121392e-7,0.869922e-11);
rd.saturated_liquid_density_in(1.8877394,0.59858531,-0.071134041,0.40327650,0.0,0.0);
rd.viscosity_liq_in(2.1818800041577616E-4,-2.36838590237531E-6,1.3282371923066061E-8,-1.0821844104196917E-10,
5.762298079646129E-13,5.502133294478409E-16,1.9437525623123426E-17,-6.687479715975144E-19,1.1172362461765993E-21);
rd.viscosity_vap_in(1.1507964713318409E-5,4.13025972870001E-8,1.998600601972917E-11,5.352228220987614E-12,
6.163698553874306E-14,-1.2615666541725562E-15,-1.3420686560142838E-17,1.3267564450833358E-19,1.2937769741708501E-21);
rd.thermal_conductivity_liq_in(0.09499678110835275,-4.7537555020566585E-4,-4.586477586766141E-7,3.241991282985186E-8,
4.850066566168331E-10,-1.1099453270713835E-11,-1.309775704321957E-13,9.740448684253897E-16,1.0670244467108473E-17);
rd.thermal_conductivity_vap_in(0.009483329607503754,5.0874888867390514E-5,-1.235230947106529E-8,2.5879027516334643E-8,
3.384514677289011E-10,-6.943747680873354E-12,-7.930356771781604E-14,6.932045271314212E-16,7.145928671486519E-18);
rd.surface_tension_in(0.011697626888583501,-1.4973507686581213E-4,1.7719089602497538E-7,5.807291774672925E-10,
3.5242874096858627E-12,-2.1565432337536383E-14,-2.595295965202523E-16,4.152703469410444E-18,3.6048844877896633E-20);
rd.saturated_vapor_density_in(-0.1125494951,
0.0,
0.0,
0.0,	
0.4892280363e3, 	
-0.1227271507e4, 	
0.1629238081e4 	,
-0.8591945110e3, 	
-0.1440371655e3, 	
0.2925609255e3,
0.0,
-0.2447993522e-2,
0.0,
0.0,	
0.1815569814e-4,
0.0,
-0.4028006124e-6, 	
0.2289562605e-7,
0.0,
0.9465371017e2,
0.0,
0.0,
0.0,
0.0);

try{gi.addElement(rd);} catch(ArrayIndexOutOfBoundsException e) 
   {System.out.println("Array Index Out Of Bounds Exception");}

/*
rd=new sogut1Data("R12","dichlorofluoromethane",
8.31442,4.616,385.2,4.1250,0.216692,384.95,120.914,339.110,-1.33018);
rd.EOS_in(2.337370817,-6.412406390,4.830185249,-1.965519089,-3.007980190,
          6.836347882,-3.379592003,3.125967270,-7.177260909,4.136982393,
          0.4704538125e-2,-0.1573667708e-1,0.1398447840,-0.1863103713,0.1221380142e-1,
          0.0,0.0,0.0,0.0,0.0,
          0.0,0.0,0.0,0.0,0.0,
          0.0,0.0,0.0);

rd.vapor_pressure_in(-6.9670659,1.6788237,-4.079537026,4.482102063,-5.064659805);
rd.ideal_gas_heat_capacity_in(0.85309e1,-0.58573e-1,0.36561e-2,-0.71739e-5,0.77377e-8,-0.35804e-11);
rd.saturated_liquid_density_in(0.3036359e1,-0.1156439e2,0.4347273e2,-0.7285118e2,0.5858101e2,-0.1794327e2);
rd.saturated_vapor_density_in(0.0,
0.0,
0.0,
0.0,
0.0,
0.7889048250e-5,
0.0,
0.0,
0.0,
0.0,
-0.1900330949,
0.0,
-0.5947813469e1,
0.0,
0.0,
0.0,
0.0,
0.0,
0.0,
0.0,
0.0,
-0.5126536438e3,
0.1142994928e4,
-0.6905751515e3);

try{gi.addElement(rd);} catch(ArrayIndexOutOfBoundsException e) 
   {System.out.println("Array Index Out Of Bounds Exception");}
*/
int i;
Enumeration n1=gi.elements();
while(n1.hasMoreElements())
{
   rd=(sogut1Data)n1.nextElement();
   s+=rd.refrigerantName+" ";
   if (rd.refrigerantName.equals(rName)) {break;}
}
return s;
}
}