def get_mf6_mshape(disfile):
    """Return the shape of the MODFLOW 6 model.

    Parameters
    ----------
    disfile : str
        path to a MODFLOW 6 discretization file

    Returns
    -------
    mshape : tuple
        tuple with the shape of the MODFLOW 6 model.

    """
    with open(disfile, "r") as f:
        lines = f.readlines()

    d = {}
    for line in lines:

        # Skip over blank and commented lines
        ll = line.strip().split()
        if len(ll) < 2:
            continue
        if line.strip()[0] in ["#", "!"]:
            continue

        for key in ["NODES", "NCPL", "NLAY", "NROW", "NCOL"]:
            if ll[0].upper() in key:
                d[key] = int(ll[1])

    if "NODES" in d:
        mshape = (d["NODES"],)
    elif "NCPL" in d:
        mshape = (d["NLAY"], d["NCPL"])
    elif "NLAY" in d:
        mshape = (d["NLAY"], d["NROW"], d["NCOL"])
    else:
        print(d)
        raise Exception("Could not determine model shape")
    return mshape