if(fvc::domainIntegrate(F[group]*rho).value() >= fuelperFG) 
{
	NewFGindex = NewFGindex + 1;
	if(group < (nf-1.5))
	{
		group++;
        	forAll(mf.boundaryFieldRef(), patchi)
        	{
        		fvPatchScalarField& ff1 = F[group-1].boundaryFieldRef()[patchi];
                	fvPatchScalarField& ff2 = F[group].boundaryFieldRef()[patchi];
                	fvPatchScalarField& mm = mf.boundaryFieldRef()[patchi];
                	forAll(mm, facei)
                	{
                	       	ff1[facei] = 0.0;
                	        ff2[facei] = mm[facei]; 
                	}
        	}
		eventRecordFile<<runTime.value()<<tab<<"New flame group is injected"<<tab<<group<<endl;		
	}
	else if(group > (nf-1.5))
	{
		//retain group, BC
		forAll(rho, celli)
		{
			F[0][celli] = F[0][celli] + F[1][celli];
			for(label k=1 ; k<nf-1; k++)
			{
				F[k][celli] = F[k+1][celli];
			}
			F[nf-1][celli] = 0.0;
		}
		
		//save current conditional profile temporarily
		scalarField QiCMC_temptemp( Ysize*(etamax+1)*nf , 0 );
		scalarField QhCMC_temptemp( (etamax+1)*nf , 0 );
		scalarField TCMC_temptemp( (etamax+1)*nf , 0);
		scalarField rhoCMC_temptemp( (etamax+1)*nf, 0);


		QiCMC_temptemp = QiCMC;
		QhCMC_temptemp = QhCMC;
		TCMC_temptemp = TCMC;
		rhoCMC_temptemp = rhoCMC;


		//set initial profile
		#include "QCMC_init.H"

		//replace conditional profile (copy conditional profile of F[k+1] to F[k])
		for(label k=0 ; k<nf-1 ; k++)
		{
			for(label j=0 ; j<=etamax ; j++)
			{
				for(label i=0 ; i<Ysize ; i++)
				{
					QiCMC[(k*(etamax+1)+j)*Ysize + i] = QiCMC_temptemp[((k+1)*(etamax+1)+j)*Ysize + i];
				}
				TCMC[j+k*(etamax+1)] = TCMC_temptemp[j+(k+1)*(etamax+1)];
				rhoCMC[j+k*(etamax+1)] = rhoCMC_temptemp[j+(k+1)*(etamax+1)];
				QhCMC[j+k*(etamax+1)] = QhCMC_temptemp[j+(k+1)*(etamax+1)];
	
			}
		}
		eventRecordFile<<runTime.value()<<tab<<"FG 0 and 1 are merged"<<tab<<"mass of FG 0 and 1 are "
		<<fvc::domainIntegrate(F[0]*rho).value()<<tab<<fvc::domainIntegrate(F[1]*rho).value()<<endl;
	}
	else
	{
		abort();
	}
}