scalarFieldArray1d C1table //AMC C1 coefficient table
(
    etamax,
    scalarField(NVar+1, 0)
);

scalarFieldArray2d Ttable
(
    etamax,
    scalarFieldArray1d
    (
        NVar+1,
        scalarField(NstList.size(), 0)
    )
);

scalarFieldArray2d Rtable
(
    etamax,
    scalarFieldArray1d
    (
        NVar+1,
        scalarField(NstList.size(), 0)
    )
);

scalarFieldArray3d Ytable
(
    Y.size(),
    scalarFieldArray2d
    (
        etamax,
        scalarFieldArray1d
        (
            NVar+1,
            scalarField(NstList.size(), 0)
        )
    )
);

Info<<"Construct C1 coefficient table"<<endl;
for(label j = 1; j < etamax - 1; j++)
{
    for(label v = 1; v < varValue.size() - 1; v++)
    {
        const scalar mf = etaValue[j];
        const scalar mfVar = varValue[v]*mf*(1.0-mf);

        CMC::BetaFunction bf(mf, mfVar);
        CMC::BetaIntegrator bi(bf, bg);

        C1table[j][v] = 1.0/bi.betaIntegrate(amcGridValue);
    }
}

Info<<"Construct Temp and Rgas table"<<endl;
for(label j=0 ; j<etamax ; j++)
{
    for(label v=0 ; v<=NVar ; v++)
    {
        const scalar mf = etaValue[j];
        const scalar mfVar = varValue[v]*mf*(1.0-mf);

        CMC::BetaFunction bf(mf, mfVar);
        CMC::BetaIntegrator bi(bf, bg);

        for(label n=0 ; n<NstList.size() ; n++)
        {
            Ttable[j][v][n] = bi.betaIntegrate(T_SLFM[n]);
            Rtable[j][v][n] = bi.betaIntegrate(Rgas_SLFM[n]);
            forAll(Y, yi)
            {
                Ytable[yi][j][v][n] = bi.betaIntegrate(Y_SLFM[yi][n]);
            }
        }
    }
    Info<<scalar(j+1)/etamax*100<<"%"<<endl;
}