multiple_pdinv does NOT return inverses

[mathDR]

Hi. In the utilities.py file, the function multiple_pdinv is described as:

def multiple_pdinv(A):
    """
    Arguments
    ---------
    A : A DxDxN numpy array (each A[:,:,i] is pd)

    Returns
    -------
    invs : the inverses of A
    hld: 0.5* the log of the determinants of A
    """

but when I call the function with random PD matrices, I don't always return the inverse.

In particular, calling the function with the PD matrix

A = np.asarray([[2.,-1.,0.],[-1.,2.,-1.],[0.,-1.,2.]])[:,:,None]

returns the correct

[[ 0.75  0.5   0.25]
 [ 0.5   1.    0.5 ]
 [ 0.25  0.5   0.75]]

but calling it with a random matrix, say

A = [[[ 1.00487128]
  [ 0.10450152]
  [ 0.78902973]]
 [[ 0.64999066]
  [ 1.92596954]
  [ 0.08218574]]
 [[ 0.59447227]
  [ 0.43842888]
  [ 1.17832435]]]

returns

invs[:,:,0] = [[ 1.4981143   0.14944808 -1.01359122]
 [-0.48097944  0.47961461  0.28862138]
 [-0.57684638 -0.2538517   1.25263636]]

but the inverse is

np.linalg.inv(A) = [[ 1.66485487 -0.40496621 -0.68925095]
 [-0.40496621  0.66577222 -0.04341129]
 [-0.68925095 -0.04341129  1.21254673]]

Am I missing an assumption on what types of tensors go into multiple_pdinv?