def METHOD_NAME():
    mu = floatX_array([-0.3, 0.7, 1.3])
    p = floatX_array([[2.0, 0, 0], [0.05, 0.1, 0], [1.0, -0.05, 5.5]])
    tau = np.dot(p, p.T)
    with pm.Model() as model:
        pm.MvNormal(
            "x",
            pt.constant(mu),
            tau=pt.constant(tau),
            initval=floatX_array([0.1, 1.0, 0.8]),
        )
    H = tau
    C = np.linalg.inv(H)
    return model.initial_point(), model, (mu, C)