Given a 3d matrix, make a block matrix, where the first and second dimensions are blocked according
    to blocksizes, and the pages are blocked using pagesizes

def get_blocks_3d(A, blocksizes, pagesizes=None):
    """
    Given a 3d matrix, make a block matrix, where the first and second dimensions are blocked according
    to blocksizes, and the pages are blocked using pagesizes
    """
    assert (A.shape[0]==A.shape[1]) and len(A.shape)==3, "can't blockify this non-square matrix, may need to use 2d version"
    N = np.sum(blocksizes)
    assert A.shape[0] == N, "bad blocksizes"
    num_blocks = len(blocksizes)
    if pagesizes == None:
        #Assume each page of A should be its own dimension
        pagesizes = range(A.shape[2])#[0]*A.shape[2]
    num_pages = len(pagesizes)
    B = np.empty(shape=(num_blocks, num_blocks, num_pages), dtype=np.object)
    count_k = 0
    #for Bk, k in enumerate(pagesizes):
    for Bk in pagesizes:
        count_i = 0
        for Bi, i in enumerate(blocksizes):
            count_j = 0
            for Bj, j in enumerate(blocksizes):
                #We want to have it count_k:count_k + k but its annoying as it makes a NxNx1 array is page sizes are set to 1
                B[Bi, Bj, Bk] = A[count_i:count_i + i, count_j:count_j + j, Bk]
                count_j += j
            count_i += i
        #count_k += k
    return B