---

import numpy as np

def broyden(f, x, eps=1e-5):
    it = 0
    e = np.eye(len(x)) * eps
    J = (f(x-e) - f(x+e)) / (-2*eps)
    x = np.array(x).reshape(len(x), 1)
    fx = f(x)

    while np.linalg.norm(fx) > eps:
        dx = -np.linalg.inv(J) @ fx
        x = x + dx
        fx = f(x)
        J = J + (fx @ dx.T) / (dx.T @ dx)
        it += 1
    return x.flatten(), fx.flatten(), it

def f(x):
    return np.array([
        1/(x[0] + 1.5) - x[1] - 1,
        (1.2 * x[0])**3 - x[1] - 1
    ])

x, y, it = broyden(f, [0.3, 0.3])
print(x, y, it, sep="\n")