def METHOD_NAME(functions, derivatives, directions, stag, msg):
    global errorlist
    for direction in directions:
        direction, fac, guards, diff_func = direction
        locations = ["centre"]
        if stag:
            locations.append(direction.lower() + "low")
        for funcs, derivative, inloc, outloc in itertools.product(
            functions, derivatives, locations, locations
        ):
            infunc, outfunc = funcs
            order, diff = derivative
            errors = []
            for nz in nzs:
                boutpp.setOption(
                    "meshD:nD".replace("D", direction),
                    "%d" % (nz + (2 * guards if direction == "x" else 0)),
                    force=True,
                )
                boutpp.setOption(
                    "meshD:dD".replace(
                        "D",
                        direction,
                    ),
                    "2*pi/(%d)" % (nz),
                    force=True,
                )
                dirnfac = direction + "*" + fac
                mesh = boutpp.Mesh(section="mesh" + direction)
                f = boutpp.create3D(infunc.replace("%s", dirnfac), mesh, outloc=inloc)
                sim = diff_func(f, method=diff, outloc=outloc)
                if sim.getLocation() != outloc:
                    cent = ["CENTRE", "CENTER"]
                    if outloc in cent and sim.getLocation() in cent:
                        pass
                    else:
                        errorlist.append(
                            "Location does not match - expected %s but got %s"
                            % (outloc, sim.getLocation())
                        )
                ana = boutpp.create3D(
                    outfunc.replace("%s", dirnfac), mesh, outloc=outloc
                )
                err = sim - ana
                err = err.getAll().flatten()
                if guards:
                    err = err[guards:-guards]
                err = np.max(np.abs(err))
                errors.append(err)
            errc = np.log(errors[-2] / errors[-1])
            difc = np.log(nzs[-1] / nzs[-2])
            conv = errc / difc
            if order - 0.1 < conv < order + 0.1:
                pass
            else:
                info = "%s - %s - %s - %s -> %s " % (
                    infunc,
                    diff,
                    direction,
                    inloc,
                    outloc,
                )
                error = "%s: %s is not working. Expected %f got %f" % (
                    msg,
                    info,
                    order,
                    conv,
                )
                errorlist.append(error)
                if doPlot:
                    from matplotlib import pyplot as plt
                    plt.plot((ana).getAll().flatten())
                    plt.plot((sim).getAll().flatten())
                    plt.show()