def _compute_axb(g1, h1, g2, h2, k):
    """
    Compute the AXB (see paper) formula
    :param np.ndarray g1: the head matrix of the first graph
    :param np.ndarray h1: the tail matrix of the first graph
    :param np.ndarray g2: the head matrix of the second graph
    :param np.ndarray h2: the tail matrix of the second graph
    :param np.ndarray k: the kernel matrix between the nodes of the first graph (lines) against the second's (columns)
    :return: the AXB product
    """
    return h1 @ (g1.T @ g2 * k) @ h2.T + h1 @ (g1.T @ h2 * k) @ g2.T +\
        g1 @ (h1.T @ g2 * k) @ h2.T + g1 @ (h1.T @ h1 * k) @ g2.T