IFstream slfmFile
		(
		        SLFMinit
		);

		word garbage, yname;
		scalar neta(0), ny(0);

	        slfmFile.getLine(garbage);//INPUT
	        slfmFile.getLine(garbage);//NO.
	        slfmFile>>neta>>ny;//number of eta and species

	        slfmFile.getLine(garbage);//blank line
	        slfmFile.getLine(garbage);//PRESSURE
	        slfmFile.getLine(garbage);//1.0
	        slfmFile.getLine(garbage);//MIXTURE

	        scalarField eta(neta, 0);
	        scalarField yi_temp(neta*(ny+1), 0), yi((etamax+1)*(ny+1),0), wi(neta*(ny+1), 0);

	        for(label j=0 ; j<neta ; j++)
	        {
	                slfmFile>>eta[j];
	        }

	        slfmFile.getLine(garbage); //blank line
	        slfmFile.getLine(garbage); //INITIAL

	        for(label i=0; i<ny ; i++) //species
	        {
	                slfmFile>>yname;
			Info<<yname<<endl;
			for(label j=0 ; j<neta ; j++)
			{
				slfmFile>>yi_temp[j*(ny+1)+i];
			}

			for(label j=0 ; j<neta ; j++)
			{
				slfmFile>>wi[j*(ny+1)+i]; //reaction rates
			}
			slfmFile.getLine(garbage);//blank line
		}
		//Temperature

		slfmFile>>yname;
		Info<<yname<<endl;
		for(label j=0 ; j<etamax+1 ; j++)
		{
			slfmFile>>yi_temp[j*(ny+1)+ny];

		}

		Interpolation(etaValue, eta, 1, ny+1, yi_temp, yi);

		for(label j = 0 ; j<=etamax ; j++)
		{
			for(label i = 0 ; i<Ysize ; i++)
			{
				QiCMC[j*Ysize+i] = yi[j*(ny+1)+i];
			}
			TCMC[j] = yi[j*(ny+1)+ny];
		}

		for(label k=1 ; k<nf ; k++)
		{
			for(label j = 0 ; j<=etamax ; j++)
			{
				for(label i = 0 ; i<Ysize ; i++)
				{
					QiCMC[(k*(etamax+1)+j)*Ysize+i] = QiCMC[j*Ysize+i];
				}
				TCMC[j+k*(etamax+1)] = TCMC[j];
			}
		}



for(label j=0 ; j <= etamax ; j++) //get conditional enthalpy 
{
        for(label k = 0 ; k<nf ; k++)
        {
                scalar Tin = TCMC[ j+k*(etamax+1) ];
                scalarField QiCMCin(Ysize, 0);

                for(label i=0 ; i<Ysize ; i++)
                {
                        QiCMCin[i] = QiCMC[(k*(etamax+1)+j)*Ysize + i];
                }

                rhoCMC[j+k*(etamax+1)] = chemistry->calculateRHOCMC(QiCMCin , Tin, Pin);

                QhCMC[j+k*(etamax+1)] = chemistry->calculateHCMC(QiCMCin , Tin, rhoCMC[ j+k*(etamax+1)], Pin);

        }
}