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

class Reaction1
{
  public int ngas;
  String reacName;
  Gas gasR[],gasP[];
  public int natom;
  Atom atomList[];
  File reacFile;
  double nr[],np[];
  double nrt,npt;

  public Reaction(String name) throws IOException
  {
  reacFile=new File("Reaction.dat");
  reacFile=new File(reacFile.getAbsolutePath());
  try
     {
     DataInputStream cfin=new DataInputStream(new FileInputStream(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 Gas[ngas];
       gasP=new Gas[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 Gas(tempGasName,nr[i]);
         gasP[i]=new Gas(tempGasName,np[i]);
         }
       }
     catch(FileNotFoundException fnfe) {System.out.println("File Not Found");}
     }   //end of constructor

     public String toString()
     {
       String s="";
       s+=reacName+"\n";
       for(int i=0;i<ngas;i++)
       {
       s+=gasR[i].toString()+" "+gasR[i].N+" "+gasP[i].N+"\n";
       }
     return s;
     }

     public String ReactiontoString()
     {
       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+=" + ";
       }
       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+=" + ";
       }
     return s;
     }


     public double H(double TR,double TP)
     {
     //reaction enthalpy
      double HH=0;
      for(int i=0;i<ngas;i++)
      {
       HH+=gasP[i].HT(TP);
       HH-=gasR[i].HT(TR);
      }
      return HH;
     }

     public double S(double TR,double TP,double pt)
     {
     //reaction entropy
      double SS=0;
      for(int i=0;i<ngas;i++)
      {
       if(np[i]!=0)
          SS+=gasP[i].S(TP,np[i]/npt*pt);
       if(nr[i]!=0)
          SS-=gasR[i].S(TR,nr[i]/nrt*pt);
      }
      return SS;
     }

     public double S(double TR,double TP)
     {
     //reaction entropy
      double SS=0;
      for(int i=0;i<ngas;i++)
      {
       if(np[i]!=0)
          SS+=gasP[i].S(TP,np[i]/npt);
       if(nr[i]!=0)
          SS-=gasR[i].S(TR,nr[i]/nrt);
      }
      return SS;
     }

     public double G(double TR,double TP,double pt)
     {
     // reaction gibbs energy
     return H(TR,TP)-TR*S(TR,TP,pt);
     }

     public double G(double TR,double TP)
     {
     // reaction gibbs energy for P=1 bar
     return H(TR,TP)-TR*S(TR,TP);
     }

     public double Taf(double TR,double xl,double 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,fxr;
     int maxit=500,iter=0;
     es=0.000001;
     ea=1.1*es;
     while((ea>es)&&(iter<maxit))
     {
     xr=(xl+xu)/2.0;
     iter++;
     if((xl+xu)!=0)
        { ea=Math.abs((xu-xl)/(xl+xu))*100;}
     fxl= H(TR,xl);
     fxr= H(TR,xr);
     test= fxl*fxr;
     if(test==0.0)       ea=0;
     else if(test<0.0)  xu=xr;
     else
       {
       xl=xr;
       }
     } //end of while
     return xr;
     }
  }


