def feed_exhaust(Uf, Ua, MS, n_CH4f, n_H2ex):
    """
    Define the settings of the model.
    The reactions are as below:
    Reforming: CH4 + H2O -> CO + 3H2
    Water gas shift: CO + H2O -> CO2 + H2
    Methane combustion: CH4 + 2O2 -> CO2 + 2H2O
    Hydrogen combustion: H2 + 1/2O2 -> H2O
    Carbon monoxide combustion: CO + 1/2O2 -> CO2
    Optional keyword arguments:
    ---------------------------
    :param Uf: Fuel ultilization, mole reductant consumed in FC per mole
               of reductant total
    :param Ua: Air ultilization, mole of air consumed in FC per mole
               of air feed
    :param MS: Methane to steam ratio, mole methane per mole water
    :param n_CH4f: mole of methane feed
    :param n_H2ex: mole of hydrogen exhaust
    :return: return the mole of methane, steam, air feed, and fuel
             ultilization
    """
    Uf = Uf  # Fuel ultilization
    Ua = Ua  # Air ultilization
    MS = MS  # Methane to steam ratio (MS): mole methane per mole water
    n_CH4f = n_CH4f
    n_H2Of = n_CH4f * MS
    n_O2f = n_CH4f * Uf * 2 / Ua
    n_N2f = n_O2f * 0.79 / 0.21
    print("mole of methane feed: " + str(n_CH4f) + " mole/s")
    print("mole of steam feed: " + str(n_H2Of) + " mole/s")
    print("mole of air feed: " + str(n_N2f + n_O2f) + " mole/s")

    n_H2ex = n_H2ex
    n_COex = n_CH4f * (1 - Uf) * 4 - n_H2ex
    n_CO2ex = n_CH4f - n_COex
    n_H2Oex = n_H2Of + 2 * n_CH4f - n_H2ex
    y_H2ex = n_H2ex / (n_H2ex + n_COex + n_CO2ex + n_H2Oex)
    y_COex = n_COex / (n_H2ex + n_COex + n_CO2ex + n_H2Oex)
    y_CO2ex = n_CO2ex / (n_H2ex + n_COex + n_CO2ex + n_H2Oex)
    y_H2Oex = n_H2Oex / (n_H2ex + n_COex + n_CO2ex + n_H2Oex)

    print("Anode exhaust: ")
    print("y_H2ex: " + str(y_H2ex))
    print("y_COex: " + str(y_COex))
    print("y_CO2ex: " + str(y_CO2ex))
    print("y_H2Oex: " + str(y_H2Oex))
    print("Total mole/s: " + str(n_H2ex + n_COex + n_CO2ex + n_H2Oex))

    n_N2ex = n_N2f
    n_O2ex = n_O2f - n_CH4f * Uf * 2
    y_O2ex = n_O2ex / (n_O2ex + n_N2ex)
    y_N2ex = n_N2ex / (n_O2ex + n_N2ex)
    print("Cathode exhaust: ")
    print("y_O2ex: " + str(y_O2ex))
    print("y_N2ex: " + str(y_N2ex))
    print("Total mole/s: " + str(n_O2ex + n_N2ex))

    return [n_CH4f, n_H2Of, n_O2f, n_N2f, Uf]