//set eta-space
scalar upperSum(0), upperDelta(0);
scalar lowerSum(0), lowerDelta(0);

//summation of geometric progression (eta > stoi)
upperSum = Foam::pow(stretchFac,upperN)-1.0;
upperSum = upperSum/(stretchFac-1.0);

//summation of geometric progression (eta < stoi)
lowerSum = Foam::pow(stretchFac,lowerN)-1.0;
lowerSum = lowerSum/(stretchFac-1.0);

//set reference delta eta
upperDelta = (1.0-stoiMF)/upperSum;
lowerDelta = stoiMF/lowerSum;

//set boundary and stoichiometric value of eta
etaValue[0] = 0.0;
etaValue[lowerN] = stoiMF;
etaValue.last() = 1.0;

//set internal value of eta (eta > stoi)
for(label j=0 ; j<upperN-1 ; j++)
{
    etaValue[lowerN+j+1] = etaValue[lowerN+j] + upperDelta*pow(stretchFac,j);
}

//set internal value of eta (eta < stoi)
for(label j=0 ; j<lowerN-1 ; j++)
{
    etaValue[lowerN-j-1] = etaValue[lowerN-j] - lowerDelta*pow(stretchFac,j);
}

Info<<etaValue<<endl;

//calculate delta eta
for(label j=1 ; j<etamax-1 ; j++)
{
    deta_h[j] = etaValue[j+1] - etaValue[j];
    deta_l[j] = etaValue[j] - etaValue[j-1];
    deta_m[j] = 0.5*(etaValue[j+1] - etaValue[j-1]);
}

//Info<<deta_h<<endl;
//Info<<deta_l<<endl;
//Info<<deta_m<<endl;

for(label j=0 ; j<etamax ; j++)
{
    etaIndex[j] = j;
}