butools.ph.CanonicalFromPH3¶

butools.ph.CanonicalFromPH3()¶

Matlab:	`[beta, B] = CanonicalFromPH3(alpha, A, prec)`
Mathematica:	`{beta, B} = CanonicalFromPH3[alpha, A, prec]`
Python/Numpy:	`beta, B = CanonicalFromPH3(alpha, A, prec)`

Returns the canonical form of an order-3 phase-type distribution.

Parameters:

alpha : matrix, shape (1,3)

Initial vector of the phase-type distribution

A : matrix, shape (3,3)

Transient generator of the phase-type distribution

prec : double, optional

Numerical precision, default value is 1e-14

Returns:

beta : matrix, shape (1,3)

The initial probability vector of the canonical form

B : matrix, shape (3,3)

Transient generator matrix of the canonical form

Notes

This procedure calculates 5 moments of the input and calls ‘PH3From5Moments’.

Examples

For Matlab:

>>> a = [0.1,0.9,0];
>>> A = [-6.2, 2, 0; 2, -9, 1; 1, 0, -3];
>>> [b, B] = CanonicalFromPH3(a, A);
>>> disp(b);
      0.58305      0.32736     0.089589
>>> disp(B);
      -9.9819            0            0
       5.3405      -5.3405            0
            0       2.8776      -2.8776
>>> Cm = SimilarityMatrix(A, B);
>>> err1 = norm(A*Cm-Cm*B);
>>> err2 = norm(a*Cm-b);
>>> disp(max(err1, err2));
   1.0269e-14

For Mathematica:

>>> a = {0.1,0.9,0};
>>> A = {{-6.2, 2, 0},{2, -9, 1},{1, 0, -3}};
>>> {b, B} = CanonicalFromPH3[a, A];
>>> Print[b];
{0.5830541253440302, 0.3273566132692404, 0.08958926138672949}
>>> Print[B];
{{-9.981920626264277, 0., 0.},
 {5.340471031780809, -5.340471031780809, 0.},
 {0., 2.8776083419564285, -2.8776083419564285}}
>>> Cm = SimilarityMatrix[A, B];
>>> err1 = Norm[A.Cm-Cm.B];
>>> err2 = Norm[a.Cm-b];
>>> Print[Max[err1, err2]];
2.7041031044230435*^-13

For Python/Numpy:

>>> a = ml.matrix([[0.1,0.9,0]])
>>> A = ml.matrix([[-6.2, 2, 0],[2, -9, 1],[1, 0, -3]])
>>> b, B = CanonicalFromPH3(a, A)
>>> print(b)
[[ 0.58305  0.32736  0.08959]]
>>> print(B)
[[-9.98192  0.       0.     ]
 [ 5.34047 -5.34047  0.     ]
 [ 0.       2.87761 -2.87761]]
>>> Cm = SimilarityMatrix(A, B)
>>> err1 = la.norm(A*Cm-Cm*B)
>>> err2 = la.norm(a*Cm-b)
>>> print(np.max(err1, err2))
8.94161968959e-14