def var_quadratic_sum(A, C, H, beta, x0):
 (A, C, H) = list(map(np.atleast_2d, (A, C, H)))
 x0 = np.atleast_1d(x0)
 Q = scipy.linalg.solve_discrete_lyapunov((sqrt(beta) * A.T), H)
 cq = dot(dot(C.T, Q), C)
 v = ((np.trace(cq) * beta) / (1 - beta))
 q0 = (dot(dot(x0.T, Q), x0) + v)
 return q0

def var_quadratic_sum ( A , C , H , beta , x0 ) : ( A , C , H ) = list ( map ( np . atleast_2d , ( A , C , H ) ) ) x0 = np . atleast_1d ( x0 ) Q = scipy . linalg . solve_discrete_lyapunov ( ( sqrt ( beta ) * A . T ) , H ) cq = dot ( dot ( C . T , Q ) , C ) v = ( ( np . trace ( cq ) * beta ) / ( 1 - beta ) ) q0 = ( dot ( dot ( x0 . T , Q ) , x0 ) + v ) return q0
computes the expected discounted quadratic sum .