'A simple algorithm for computing the robust policy F and the
corresponding value function P, based around straightforward
iteration with the robust Bellman operator.  This function is
easier to understand but one or two orders of magnitude slower
than self.robust_rule().  For more information see the docstring
of that method.
Parameters
P_init : array_like(float, ndim=2), optional(default=None)
The initial guess for the value function matrix.  It will
be a matrix of zeros if no guess is given
max_iter : scalar(int), optional(default=80)
The maximum number of iterations that are allowed
tol : scalar(float), optional(default=1e-8)
The tolerance for convergence
Returns
F : array_like(float, ndim=2)
The optimal control matrix from above
P : array_like(float, ndim=2)
The positive semi-definite matrix defining the value
function
K : array_like(float, ndim=2)
the worst-case shock matrix K, where
:math:`w_{t+1} = K x_t` is the worst case shock'

def robust_rule_simple(self, P_init=None, max_iter=80, tol=1e-08):

    (A, B, C, Q, R) = (self.A, self.B, self.C, self.Q, self.R)
    (beta, theta) = (self.beta, self.theta)
    P = (np.zeros((self.n, self.n)) if (P_init is None) else P_init)
    (iterate, e) = (0, (tol + 1))
    while ((iterate < max_iter) and (e > tol)):
        (F, new_P) = self.b_operator(self.d_operator(P))
        e = np.sqrt(np.sum(((new_P - P) ** 2)))
        iterate += 1
        P = new_P
    I = np.identity(self.j)
    S1 = P.dot(C)
    S2 = C.T.dot(S1)
    K = inv(((theta * I) - S2)).dot(S1.T).dot((A - B.dot(F)))
    return (F, K, P)
