def METHOD_NAME(phase, f, **kwargs):
    r"""
    Wrapper function for calling models in the ``chemicals`` package
    Parameters
    ----------
    phase : dict
        The OpenPNM Species object for which this model should calculate
        values. This object should ideally have all the necessary chemical
        properties in its ``params`` attribute, although this is optional as
        discussed below in relation to ``kwargs``.
    f : function
        The handle of the function to apply, such as
        ``chemicals.viscosity.Letsou_Stiel``.
    kwargs
        By default this function will use
        ``openpnm.models.phase.default_argmap`` to determine which names on
        ``phase`` correspond to each argument required by ``f``. For
        instance, ``default_argmap['Tc'] = 'param.critical_temperature'``
        so any functions that need ``Tc`` will receive
        ``phase['param.critical_temperature']``. Some models only require
        the normal thermodynamic parameters, like ``T``, ``Tc``, and ``Pc``,
        so in many cases no additional arguments are needed beyond whats in the
        default argument map. However, other models require values that are
        themselves calcualted by another model, like ``Cvm``, or perhaps a
        list of constants like ``a0``, ``a1``, etc. It is necessary to provide
        these as keyword arguments, and they will be included in
        ``default_argmap`` as new entries or overwriting existing ones.
        So ``mu='pore.blah'`` will pass ``phase['pore.blah']`` to the
        ``mu`` argument of ``f``. Any arguments ending with an ``s`` are
        assumed to refer to the individual properties of pure components in a
        mixture. So ``mus`` means fetch the viscosity of each component as a
        list, like ``[1.3e-5, 1.9e-5]``. This function will trim the trailing
        ``s`` off any argument name before checking the ``default_argmap`` so
        that the normal pure component values can be looked up.
    Notes
    -----
    This wrapper works with both pure and mixture phases, but the mixture
    models are very slow due to the way ``chemicals`` vectorizes code. For pure
    species it allows the computation of values at many different conditions
    in a vectorized way. The means that the conditions in each pore, such as
    temperature, pressure, etc can be passed and iterpreted as a list of
    conditions. For mixture models, however, the vectorization is done over
    the compositions, at a *fixed* condition. This means that we must do a
    pure-python for-loop (ie. slow) for each individual pore. As such, we
    have re-implemented several of the most useful mixing models offered by
    ``chemicals`` in ``OpenPNM``, and include unit tests to ensure both
    implementations agree.
    """
    import chemicals as _chemicals
    # Update default argmap with any user supplied values
    argmap = default_argmap.copy()
    for k, v in kwargs.items():
        argmap[k] = v
    args = _get_items_from_target(phase, f, argmap)
    # f = getattr(_chemicals.numba_vectorized, f.__name__)
    msg = f"Numba version failed for {f.__name__}, reverting to pure python"
    if len(set(['xs', 'yz', 'zs']).intersection(args.keys())):
        # Call function in for-loop for each pore since they are not vectorized
        logger.info(msg)
        vals = _np.zeros(phase.Np)
        for pore in phase.Ps:
            a = {}
            for item in args.keys():
                if item.endswith('s'):
                    try:
                        a[item] = [args[item][i][pore] for i in range(len(args[item]))]
                    except TypeError:
                        a[item] = args[item]
                else:
                    a[item] = args[item][pore]
            vals[pore] = f(*list(a.values()))
    else:
        try:
            # Get the numba vectorized version of f, or else numpy arrays don't work
            f = getattr(_chemicals.numba_vectorized, f.__name__)
            vals = f(*args.values())
        except AssertionError:
            f = getattr(_chemicals, f.__name__)
            vals = f(*args.values())
            logger.warn(msg)
    return vals