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

class Reaction
{
  //======================================================//
  // Dr. Turhan Coban                                     //
  // EGE University School of Engineering,                //
  // Mechanical Engineering Department                    //
  // Phone : 90 (232) 3434000-5387                        //
  // email : turhan.coban@ege.edu.tr                      //
  //======================================================//
  // Class Description :
  //------------------------------------------------------
  // This program calculates energy, adiabatic flame temperature
  // entropy and gibbs free energy for a given reaction
  // Basic reaction data is given in the file called Reaction.dat
  // Format of the data :
  // reactionname
  // number_of_the_chemicals
  // gas_name reactant_mole prodoct_mole
  // .......
  // gas_name reactant_mole prodoct_mole
  //-----------------
  // sample data :
  //-----------------
  // c2h6combustion
  // 5
  // c2h6 1 0
  // o2 3 0
  // n2 11.28 11.28
  // co2 0 2
  // h2o 0 2
  //------------------
  // the available gas_name's can be checked from data file
  // Gas.dat
  // Variable Identification
  //-----------------------------------------------------
  // ngas     : number of gases in reaction (integer variable)
  // reacname : reaction name  (String variable)
  // gasR,gasP: gas objects for reactants and products
  //            (array objects of class Gas)
  //            see Gas.java for details
  // natom    : number of atoms in the reaction
  // nr,np    : number of Moles of reactants and products
  //            (double vector array)
  // nrt,npt  : total number of moles for reactants and products
  // reacFile : reaction file channel variable (Class File)
  // atomList : list of the atoms in the reaction
  public boolean mole; //switch for mole/mass basis
  public boolean SI;   //switch for SI/EN units
  public int ngas;
  String reacName;
  Gas1 gasR[],gasP[];
  public int natom;
  Atom atomList[];
  File reacFile;
  double nr[],np[];
  double nrt,npt;

  public Reaction(String name) throws IOException
  {
  mole=true;
  SI=true;
  reacFile=new File("Reaction.txt");
  reacFile=new File(reacFile.getAbsolutePath());
  try
     {    
     //DataInputStream cfin=new DataInputStream(new FileInputStream(reacFile));
     BufferedReader cfin=new BufferedReader(new FileReader(reacFile));
     int ierror=1;
     natom=0;
     try
       {
       while(cfin!=null)
         {
         reacName=Text.readString(cfin);
         if(reacName.equals(name)) {ierror=0;break;}
         }
       }
       catch(EOFException e_eof)
         {
         System.out.println("error required gas mixture "+name+" is not found");
         cfin.close();return;
         }
       ngas=Text.readInt(cfin);
       gasR=new Gas1[ngas];
       gasP=new Gas1[ngas];
       nr=new double[ngas];
       np=new double[ngas];
       nrt=0;
       npt=0;
       String tempGasName;
       Double n1,n2;
       for(int i=0;i<ngas;i++)
         {
         tempGasName=Text.readString(cfin);
         nr[i]=Text.readDouble(cfin);
         np[i]=Text.readDouble(cfin);
         nrt+=nr[i];
         npt+=np[i];
         gasR[i]=new Gas1(tempGasName,nr[i]);
         gasP[i]=new Gas1(tempGasName,np[i]);
         }
       }
     catch(FileNotFoundException fnfe) {System.out.println("File Not Found");}
     arrange_atoms();
     }   //end of constructor

     public String readReactionNames() throws IOException
    {
    String temp1="";
    reacFile=new File("Reaction.txt");
    reacFile=new File(reacFile.getAbsolutePath());
    try
    {
    //DataInputStream cfin=new DataInputStream(new FileInputStream(reacFile));
    BufferedReader cfin=new BufferedReader(new FileReader(reacFile));
    int ierror=1;
    natom=0;
    try
       {
       while(cfin!=null)
         {
         reacName=Text.readString(cfin);
         temp1=temp1+reacName+ " ";
         ngas=Text.readInt(cfin);
         double nri;
         double npi;
         String tempGasName;
         Double n1,n2;
            for(int i=0;i<ngas;i++)
            {
            tempGasName=Text.readString(cfin);
            nri=Text.readDouble(cfin);
            npi=Text.readDouble(cfin);
            }  //end try
         }   //end while
       }  catch(EOFException e_eof)  {cfin.close();} // end of first try
     } // end of second try
     catch(FileNotFoundException fnfe) {System.out.println("File Not Found");}
     return temp1;
}


     public void arrange_atoms() throws IOException
     {
     // atomic structure of the reactants
     int i,j;
        for(i=0;i<ngas;i++)
        {
             for(j=0;j<gasR[i].natom;j++)
	     {
             add_atom(i,j);
             }
        }
     }

     public int  add_atom(int i,int j) throws IOException
     {
     //determine and order the atomic structure
     int k;
       for(k=0;k<natom;k++)
       {
       if(gasR[i].atomList[j].symbol.equals(atomList[k].symbol))
         {
         atomList[k]=new Atom(atomList[k].symbol,atomList[k].N+gasR[i].atomList[j].N*gasR[i].N);
         return 1;
         }
       }
     Atom atomL[];
     atomL=new Atom[natom+1];
     for(k=0;k<natom;k++)
       atomL[k]=new Atom(atomList[k]);
     atomL[natom]=new Atom(gasR[i].atomList[j].symbol,gasR[i].atomList[j].N*gasR[i].N);
     atomList=atomL;
     natom+=1;
     return 2;
     }

     public String toString()
     {
       // writes chemical reaction in chemistry norm
       String s="";
       for(int i=0;i<ngas;i++)
       {
       if(nr[i]==1)
           { s+=gasR[i].toString(); }
       else if(nr[i]!=0)
           { s+=nr[i]+" "+gasR[i].toString(); }
       if(i!=(ngas-1) && nr[i]!=0 ) s+=" + ";
       }
       if(s.endsWith(" + "))
         s=s.substring(0,s.length()-3);
       s+=" --> ";
       for(int i=0;i<ngas;i++)
       {
       if(np[i]==1)
           { s+=gasP[i].toString(); }
       else if(np[i]!=0)
           { s+=np[i]+" "+gasP[i].toString(); }
       if(i!=(ngas-1) && np[i]!=0 ) s+=" + ";
       }
       if(s.endsWith(" + "))
         s=s.substring(0,s.length()-3);
     return s;
     }

     public String toString(String ch)
     {
     //return the c
     String s="";
     int i,j;
     if(ch.equals("name"))
        s=s+reacName+"\n";
     else if(ch.equals("formula"))
        {
        s=toString();
        }
     else if(ch.equals("composition"))
        {
        for(i=0;i<natom;i++)
        s=s+atomList[i].toString()+" ";
        }
     return s;
     }

     public double H(double TR,double TP)
     {
     double TP1,TR1;
     if(!SI) {TR1=TR/1.8;TP1=TP/1.8;}
     else {TR1=TR;TP1=TP;}
     // reaction enthalpy
      double HH=0;
      double h1,h2;
      for(int i=0;i<ngas;i++)
      {
       h1=gasP[i].HT(TP1);
       h2=gasR[i].HT(TR1);
	   HH+=h1;
       HH-=h2;
       //System.out.println(gasP[i].gasName+"product  h1= "+h1);
       //System.out.println(gasR[i].gasName+"reactant "+h2);
      }
      if(SI)
        return HH;
      else
        return HH/0.42992;
     }

     public double S(double TR,double TP,double pt)
     {
     //reaction entropy
      double TP1,TR1,pt1;
      double s1,s2;
      if(!SI) {TR1=TR/1.8;TP1=TP/1.8;pt1=pt/14.504;}
      else {TR1=TR;TP1=TP;pt1=pt;}
      double SS=0;
      for(int i=0;i<ngas;i++)
      {
       if(np[i]!=0)
          SS+=gasP[i].Si(TP1,np[i]/npt*pt1);
       if(nr[i]!=0)
          SS-=gasR[i].Si(TR1,nr[i]/nrt*pt1);
      }
      if(SI)
        return SS;
      else
        return SS/0.238846;
     }

     public double S(double TR,double TP)
     {
     //reaction entropy
      double TP1,TR1,pt1,pt;
      if(!SI) {pt=14.696;TR1=TR/1.8;TP1=TP/1.8;pt1=pt/14.504;}
      else {pt=1.0;TR1=TR;TP1=TP;pt1=pt;}
      double SS=0;
      for(int i=0;i<ngas;i++)
      {
       if(np[i]!=0)
          SS+=gasP[i].Si(TP1,np[i]/npt*pt1);
       if(nr[i]!=0)
          SS-=gasR[i].Si(TR1,nr[i]/nrt*pt1);
      }
      if(SI)
        return SS;
      else
        return SS/0.238846;
     }

     public double G(double TR,double TP,double pt)
     {
     // reaction gibbs energy
      double TP1,TR1,pt1;
      if(!SI) {TR1=TR/1.8;TP1=TP/1.8;pt1=pt/14.504;}
      else {TR1=TR;TP1=TP;pt1=pt;}
      double GG=0;
      for(int i=0;i<ngas;i++)
      {
       if(np[i]!=0)
          GG+=gasP[i].GT(TP1,np[i]/npt*pt1);
       if(nr[i]!=0)
          GG-=gasR[i].GT(TR1,nr[i]/nrt*pt1);
      }
      if(SI)
        return GG;
      else
        return GG/0.42992;
     }

     public double G(double TR,double TP)
     {
     // reaction gibbs energy
      double TP1,TR1,pt1,pt;
      if(!SI) {pt=14.696;TR1=TR/1.8;TP1=TP/1.8;pt1=pt/14.504;}
      else {pt=1.0;TR1=TR;TP1=TP;pt1=pt;}
      double GG=0;
      for(int i=0;i<ngas;i++)
      {
       if(np[i]!=0)
          GG+=gasP[i].GT(TP1,np[i]/npt*pt1);
       if(nr[i]!=0)
          GG-=gasR[i].GT(TR1,nr[i]/nrt*pt1);
      }
      if(SI)
        return GG;
      else
        return GG/0.42992;
     }

     public double G0(double TR,double TP)
     {
     // reaction gibbs energy
      double TP1,TR1;
      if(!SI) {TR1=TR/1.8;TP1=TP/1.8;}
      else {TR1=TR;TP1=TP;}
      double GG=0;
      for(int i=0;i<ngas;i++)
      {
       if(np[i]!=0)
          GG+=gasP[i].GT(TP1);
       if(nr[i]!=0)
          GG-=gasR[i].GT(TR1);
      }
      if(SI)
        return GG;
      else
        return GG/0.42992;
     }

     public double K(double TR,double TP)
     {
     // Equilibrium costant. Usually TR and TP will be in
     // the same temperature
      double TP1,TR1;
      if(!SI) {TR1=TR/1.8;TP1=TP/1.8;}
      else {TR1=TR;TP1=TP;}
     return Math.exp(G0(TR1,TP1)/(-8.3145*TP1));
     }

     public double Taf(double TR,double xl,double xu)
     {
      double TR1,xl1,xu1;
      if(!SI) {TR1=TR/1.8;xl1=xl/1.8;xu1=xu/1.8;}
      else {TR1=TR;xl1=xl;xu1=xu;}
     // Adyabatic flame temperature using Bisection method
     // bisection method to find roots of zero energy change
     // Q=0
     // defination of variables :
     // xl : lower guess
     // xu : upper guess
     // xr : root estimate
     // es : stopping criterion
     // ea :approximate error
     // maxit : maximum iterations
     // iter  : number of iteration
     double test;
     double xr=0;
     double es,ea;
     double fxl=0,fxr=0;
     int maxit=500,iter=0;
     es=0.000001;
     ea=1.1*es;
     while((ea>es)&&(iter<maxit))
     {
     xr=(xl1+xu1)/2.0;
     iter++;
     if((xl1+xu1)!=0)
        { ea=Math.abs((xu1-xl1)/(xl1+xu1))*100;}
     fxl= H(TR,xl1);
     fxr= H(TR,xr);
     test= fxl*fxr;
     if(test==0.0)       ea=0;
     else if(test<0.0)  xu1=xr;
     else
       {
       xl1=xr;
       }
     } //end of while
      if(SI)
        return xr;
      else
        return xr*1.8;
     }
     
     double Vi(double TR,double TP,int Ni)
     {
     return -G(TR,TP)/(96487*Ni);
	 }
	 
	 double Vi0(double TR,double TP,int Ni)
     {
     return -G0(TR,TP)/(96487*Ni);
	 }

public boolean base(String s)
{
if(s.equals("mole")) {mole=true;}
else                 {mole=false;}
for(int i=0;i<ngas;i++)
{
gasR[i].base("mole");
gasP[i].base("mole");
}
return mole;
}

public boolean unit(String s)
{
if(s.equals("SI")) {SI=true;}
else               {SI=false;}
for(int i=0;i<ngas;i++)
{
gasR[i].unit("SI");
gasP[i].unit("SI");
}
return SI;
}



}


