scalarField AMCcoeff(CMC::AMC(etaValue));

//cell loop

forAll(mf, cellI)
{


    scalar jl(0), jh(0), vl(0), vh(0);
    scalar jfac(0), vfac(0);
    scalar C1coeff(0);

    //find eta-index, factor
    if(mf[cellI] < etaValue[1])
    {
        jl = 0;
        jh = 1;
    }
    else if(mf[cellI] > etaValue[etamax-2])
    {
        jl = etamax-2;
        jh = etamax-1;
    }
    else
    {
        jl = label( interpolateXY(mf[cellI], etaValue, etaIndex) );
        jh = jl+1;
    }
    jfac = max(0.0, (mf[cellI]-etaValue[jl])/(etaValue[jh]-etaValue[jl]));

    //find var-index, factor
    scalar scaledVar = min(0.99999, mfVar[cellI]/(mf[cellI]*(1.0-mf[cellI])+SMALL));
    if(scaledVar < varValue[1])
    {
        vl = 0;
        vh = 1;
    }
    else if(scaledVar > varValue[NVar-1])
    {
        vl = NVar-1;
        vh = NVar;
    }
    else
    {
        vl = label( interpolateXY(scaledVar, varValue, varIndex) );
        vh = vl+1;
    }
    vfac = max(0.0, (scaledVar-varValue[vl])/(varValue[vh]-varValue[vl]));


    //Bi-linear interpolation on j and v
    //Numerical recipes, 2nd Ed. p.117
    C1coeff
    = C1table[jl][vl]*(1-jfac)*(1-vfac)
    + C1table[jh][vl]*(jfac)*(1-vfac)
    + C1table[jh][vh]*(jfac)*(vfac)
    + C1table[jl][vh]*(1-jfac)*(vfac);

    jlc[cellI] = jl;
    jhc[cellI] = jh;
    jfc[cellI] = jfac;
    vlc[cellI] = vl;
    vhc[cellI] = vh;
    vfc[cellI] = vfac;


    scalar Coeff = AMCcoeff[lowerN];

    if (Equilibrium == true)
    {
        Neta[lowerN][cellI] = 0; //Equilibrium
    }
    else
    {
        Neta[lowerN][cellI] = mfVar[cellI]*epsk[cellI]*Coeff*C1coeff; //CSDR at stoichiometic m.f.
    }
}