butools.ph.MinimalRepFromME¶
-
butools.ph.MinimalRepFromME()¶ Matlab: [beta, B] = MinimalRepFromME(alpha, A, how, precision)Mathematica: {beta, B} = MinimalRepFromME[alpha, A, how, precision]Python/Numpy: beta, B = MinimalRepFromME(alpha, A, how, precision)Returns the minimal representation of the given ME distribution.
Parameters: alpha : vector, shape (1,M)
The initial vector of the matrix-exponential distribution.
A : matrix, shape (M,M)
The matrix parameter of the matrix-exponential distribution.
how : {“obs”, “cont”, “obscont”, “moment”}, optional
Determines how the representation is minimized. Possibilities: ‘obs’: observability, ‘cont’: controllability, ‘obscont’: the minimum of observability and controllability order, ‘moment’: moment order (which is the default).
precision : double, optional
Precision used by the Staircase algorithm. The default value is 1e-12.
Returns: beta : vector, shape (1,N)
The initial vector of the minimal representation
B : matrix, shape (N,N)
The matrix parameter of the minimal representation
References
[R34] P. Buchholz, M. Telek, “On minimal representation of rational arrival processes.” Madrid Conference on Qeueuing theory (MCQT), June 2010. Examples
For Matlab:
>>> a = [1.0/6,1.0/6,1.0/6,1.0/6,1.0/6,1.0/6]; >>> A = [-1., 0., 0., 0., 0., 0.; 0.5, -2., 1., 0., 0., 0.; 1., 0., -3., 1., 0., 0.; 1., 0., 1., -4., 1., 0.; 4., 0., 0., 0., -5., 0.; 5., 0., 0., 0., 0., -6.]; >>> [b, B] = MinimalRepFromME(a, A, 'cont'); >>> disp(b); 1 1.3878e-16 >>> disp(B); -1.4011 0.48448 0.49585 -1.5989 >>> [b, B] = MinimalRepFromME(a, A, 'obs'); >>> disp(b); 0.16667 0.16667 0.16667 0.16667 0.16667 0.16667 >>> disp(B); -1 0 0 0 0 0 0.5 -2 1 0 0 0 1 0 -3 1 0 0 1 0 1 -4 1 0 4 0 0 0 -5 0 5 0 0 0 0 -6 >>> [b, B] = MinimalRepFromME(a, A, 'obscont'); >>> disp(b); 1 1.3878e-16 >>> disp(B); -1.4011 0.48448 0.49585 -1.5989 >>> [b, B] = MinimalRepFromME(a, A, 'moment'); >>> disp(b); 0.5 0.5 >>> disp(B); -2.52 1.6467 -0.48 -0.48 >>> a = [2.0/3,1.0/3]; >>> A = [-1., 1.; 0., -3.]; >>> [b, B] = MinimalRepFromME(a, A, 'cont'); >>> disp(b); 0.66667 0.33333 >>> disp(B); -1 1 0 -3 >>> [b, B] = MinimalRepFromME(a, A, 'obs'); >>> disp(b); 1 >>> disp(B); -1 >>> [b, B] = MinimalRepFromME(a, A, 'obscont'); >>> disp(b); 1 >>> disp(B); -1 >>> [b, B] = MinimalRepFromME(a, A, 'moment'); >>> disp(b); 1 >>> disp(B); -1 >>> b = [0.2,0.3,0.5]; >>> B = [-1., 0., 0.; 0., -3., 1.; 0., -1., -3.]; >>> [a, A] = MonocyclicPHFromME(b, B); >>> disp(a); Columns 1 through 6 0.0055089 0.0090301 0.016938 0.015216 0.0053543 0.0087356 Columns 7 through 9 0.052486 0.22657 0.66016 >>> disp(A); Columns 1 through 6 -1 1 0 0 0 0 0 -2.4226 2.4226 0 0 0 0 0 -2.4226 2.4226 0 0 0 0.26232 0 -2.4226 2.1603 0 0 0 0 0 -4.2414 4.2414 0 0 0 0 0 -4.2414 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 Columns 7 through 9 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 4.2414 0 0 -4.2414 4.2414 0 0 -4.2414 4.2414 0 0 -4.2414 >>> [b, B] = MinimalRepFromME(a, A, 'cont'); >>> disp(b); Columns 1 through 6 0.0055089 0.0090301 0.016938 0.015216 0.0053543 0.0087356 Columns 7 through 9 0.052486 0.22657 0.66016 >>> disp(B); Columns 1 through 6 -1 1 0 0 0 0 0 -2.4226 2.4226 0 0 0 0 0 -2.4226 2.4226 0 0 0 0.26232 0 -2.4226 2.1603 0 0 0 0 0 -4.2414 4.2414 0 0 0 0 0 -4.2414 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 Columns 7 through 9 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 4.2414 0 0 -4.2414 4.2414 0 0 -4.2414 4.2414 0 0 -4.2414 >>> [b, B] = MinimalRepFromME(a, A, 'obs'); >>> disp(b); 1 2.0817e-17 -5.5511e-17 >>> disp(B); -2.8362 0.036222 -4.4409e-16 -16.61 -3.3369 16.042 1.1643 -0.051724 -0.82688 >>> Cm = SimilarityMatrix(B, A); >>> err1 = norm(B*Cm-Cm*A); >>> err2 = norm(b*Cm-a); >>> disp(max(err1, err2)); 9.334e-15 >>> [b, B] = MinimalRepFromME(a, A, 'obscont'); >>> disp(b); 1 2.0817e-17 -5.5511e-17 >>> disp(B); -2.8362 0.036222 -4.4409e-16 -16.61 -3.3369 16.042 1.1643 -0.051724 -0.82688 >>> Cm = SimilarityMatrix(B, A); >>> err1 = norm(B*Cm-Cm*A); >>> err2 = norm(b*Cm-a); >>> disp(max(err1, err2)); 9.334e-15 >>> [b, B] = MinimalRepFromME(a, A, 'moment'); >>> disp(b); 0.33333 0.33333 0.33333 >>> disp(B); -2.1905 1.9222 -3.3698 -1.0769 -2.3906 0.83162 -0.51037 0.8033 -2.4189 >>> Cm = SimilarityMatrix(B, A); >>> err1 = norm(B*Cm-Cm*A); >>> err2 = norm(b*Cm-a); >>> disp(max(err1, err2)); 5.5343e-15
For Mathematica:
>>> a = {1.0/6,1.0/6,1.0/6,1.0/6,1.0/6,1.0/6}; >>> A = {{-1., 0., 0., 0., 0., 0.},{0.5, -2., 1., 0., 0., 0.},{1., 0., -3., 1., 0., 0.},{1., 0., 1., -4., 1., 0.},{4., 0., 0., 0., -5., 0.},{5., 0., 0., 0., 0., -6.}}; >>> {b, B} = MinimalRepFromME[a, A, "cont"]; >>> Print[b]; {0.9999999999999999, -1.6653345369377348*^-16} >>> Print[B]; {{-1.4011480617916212, 0.48448139512495414}, {0.49584627341672904, -1.5988519382083795}} >>> {b, B} = MinimalRepFromME[a, A, "obs"]; >>> Print[b]; {0.16666666666666666, 0.16666666666666666, 0.16666666666666666, 0.16666666666666666, 0.16666666666666666, 0.16666666666666666} >>> Print[B]; {{-1., 0., 0., 0., 0., 0.}, {0.5, -2., 1., 0., 0., 0.}, {1., 0., -3., 1., 0., 0.}, {1., 0., 1., -4., 1., 0.}, {4., 0., 0., 0., -5., 0.}, {5., 0., 0., 0., 0., -6.}} >>> {b, B} = MinimalRepFromME[a, A, "obscont"]; >>> Print[b]; {0.9999999999999999, -1.6653345369377348*^-16} >>> Print[B]; {{-1.4011480617916212, 0.48448139512495414}, {0.49584627341672904, -1.5988519382083795}} >>> {b, B} = MinimalRepFromME[a, A, "moment"]; >>> Print[b]; {1/2, 1/2} >>> Print[B]; {{-2.5200000000000364, 1.6466666666667047}, {-0.48, -0.48000000000000015}} >>> a = {2.0/3,1.0/3}; >>> A = {{-1., 1.},{0., -3.}}; >>> {b, B} = MinimalRepFromME[a, A, "cont"]; >>> Print[b]; {0.6666666666666666, 0.3333333333333333} >>> Print[B]; {{-1., 1.}, {0., -3.}} >>> {b, B} = MinimalRepFromME[a, A, "obs"]; >>> Print[b]; {0.9999999999999999} >>> Print[B]; {{-1.}} >>> {b, B} = MinimalRepFromME[a, A, "obscont"]; >>> Print[b]; {0.9999999999999999} >>> Print[B]; {{-1.}} >>> {b, B} = MinimalRepFromME[a, A, "moment"]; >>> Print[b]; {1} >>> Print[B]; {{-1.}} >>> b = {0.2,0.3,0.5}; >>> B = {{-1., 0., 0.},{0., -3., 1.},{0., -1., -3.}}; >>> {a, A} = MonocyclicPHFromME[b, B]; >>> Print[a]; {0.00550893408977846, 0.00903007832853331, 0.016937512518639578, 0.015215980106503445, 0.005354337535618665, 0.008735592607040744, 0.05248568615571608, 0.22657249403204927, 0.6601593846261203} >>> Print[A]; {{-1., 1., 0., 0., 0., 0., 0., 0., 0.}, {0., -2.4226497308103743, 2.4226497308103743, 0., 0., 0., 0., 0., 0.}, {0., 0., -2.4226497308103743, 2.4226497308103743, 0., 0., 0., 0., 0.}, {0., 0.2623172489622428, 0., -2.4226497308103743, 2.1603324818481315, 0., 0., 0., 0.}, {0., 0., 0., 0., -4.241399978863847, 4.241399978863847, 0., 0., 0.}, {0., 0., 0., 0., 0., -4.241399978863847, 4.241399978863847, 0., 0.}, {0., 0., 0., 0., 0., 0., -4.241399978863847, 4.241399978863847, 0.}, {0., 0., 0., 0., 0., 0., 0., -4.241399978863847, 4.241399978863847}, {0., 0., 0., 0., 0., 0., 0., 0., -4.241399978863847}} >>> {b, B} = MinimalRepFromME[a, A, "cont"]; >>> Print[b]; {0.00550893408977846, 0.00903007832853331, 0.016937512518639578, 0.015215980106503445, 0.005354337535618665, 0.008735592607040744, 0.05248568615571608, 0.22657249403204927, 0.6601593846261203} >>> Print[B]; {{-1., 1., 0., 0., 0., 0., 0., 0., 0.}, {0., -2.4226497308103743, 2.4226497308103743, 0., 0., 0., 0., 0., 0.}, {0., 0., -2.4226497308103743, 2.4226497308103743, 0., 0., 0., 0., 0.}, {0., 0.2623172489622428, 0., -2.4226497308103743, 2.1603324818481315, 0., 0., 0., 0.}, {0., 0., 0., 0., -4.241399978863847, 4.241399978863847, 0., 0., 0.}, {0., 0., 0., 0., 0., -4.241399978863847, 4.241399978863847, 0., 0.}, {0., 0., 0., 0., 0., 0., -4.241399978863847, 4.241399978863847, 0.}, {0., 0., 0., 0., 0., 0., 0., -4.241399978863847, 4.241399978863847}, {0., 0., 0., 0., 0., 0., 0., 0., -4.241399978863847}} >>> {b, B} = MinimalRepFromME[a, A, "obs"]; >>> Print[b]; {0.9999999999999998, -6.938893903907228*^-18, 5.551115123125783*^-17} >>> Print[B]; {{-2.8362221259944986, 0.03622212599450128, 3.3306690738754696*^-16}, {-16.609614674774082, -3.3369017729049943, 16.04221903548223}, {1.1643450726193003, -0.05172359141706952, -0.8268761011005079}} >>> Cm = SimilarityMatrix[B, A]; >>> err1 = Norm[B.Cm-Cm.A]; >>> err2 = Norm[b.Cm-a]; >>> Print[Max[err1, err2]]; 1.805924023501617*^-14 >>> {b, B} = MinimalRepFromME[a, A, "obscont"]; >>> Print[b]; {0.9999999999999998, -6.938893903907228*^-18, 5.551115123125783*^-17} >>> Print[B]; {{-2.8362221259944986, 0.03622212599450128, 3.3306690738754696*^-16}, {-16.609614674774082, -3.3369017729049943, 16.04221903548223}, {1.1643450726193003, -0.05172359141706952, -0.8268761011005079}} >>> Cm = SimilarityMatrix[B, A]; >>> err1 = Norm[B.Cm-Cm.A]; >>> err2 = Norm[b.Cm-a]; >>> Print[Max[err1, err2]]; 1.805924023501617*^-14 >>> {b, B} = MinimalRepFromME[a, A, "moment"]; >>> Print[b]; {1/3, 1/3, 1/3} >>> Print[B]; {{-2.1904761904762835, 1.9221904761901123, -3.369809523809149}, {-1.0769308703296558, -2.390597900927728, 0.8316243212940291}, {-0.5103707169719298, 0.8032963136261424, -2.418925908595615}} >>> Cm = SimilarityMatrix[B, A]; >>> err1 = Norm[B.Cm-Cm.A]; >>> err2 = Norm[b.Cm-a]; >>> Print[Max[err1, err2]]; 1.6327303174900902*^-14
For Python/Numpy:
>>> a = ml.matrix([[1.0/6,1.0/6,1.0/6,1.0/6,1.0/6,1.0/6]]) >>> A = ml.matrix([[-1., 0., 0., 0., 0., 0.],[0.5, -2., 1., 0., 0., 0.],[1., 0., -3., 1., 0., 0.],[1., 0., 1., -4., 1., 0.],[4., 0., 0., 0., -5., 0.],[5., 0., 0., 0., 0., -6.]]) >>> b, B = MinimalRepFromME(a, A, "cont") >>> print(b) [[ 1.00000e+00 2.08167e-16]] >>> print(B) [[-1.40115 0.48448] [ 0.49585 -1.59885]] >>> b, B = MinimalRepFromME(a, A, "obs") >>> print(b) [[ 0.16667 0.16667 0.16667 0.16667 0.16667 0.16667]] >>> print(B) [[-1. 0. 0. 0. 0. 0. ] [ 0.5 -2. 1. 0. 0. 0. ] [ 1. 0. -3. 1. 0. 0. ] [ 1. 0. 1. -4. 1. 0. ] [ 4. 0. 0. 0. -5. 0. ] [ 5. 0. 0. 0. 0. -6. ]] >>> b, B = MinimalRepFromME(a, A, "obscont") >>> print(b) [[ 1.00000e+00 2.08167e-16]] >>> print(B) [[-1.40115 0.48448] [ 0.49585 -1.59885]] >>> b, B = MinimalRepFromME(a, A, "moment") >>> print(b) [[ 0.5 0.5]] >>> print(B) [[-2.52 1.64667] [-0.48 -0.48 ]] >>> a = ml.matrix([[2.0/3,1.0/3]]) >>> A = ml.matrix([[-1., 1.],[0., -3.]]) >>> b, B = MinimalRepFromME(a, A, "cont") >>> print(b) [[ 0.66667 0.33333]] >>> print(B) [[-1. 1.] [ 0. -3.]] >>> b, B = MinimalRepFromME(a, A, "obs") >>> print(b) [[ 1.]] >>> print(B) [[-1.]] >>> b, B = MinimalRepFromME(a, A, "obscont") >>> print(b) [[ 1.]] >>> print(B) [[-1.]] >>> b, B = MinimalRepFromME(a, A, "moment") >>> print(b) [[ 1.]] >>> print(B) [[-1.]] >>> b = ml.matrix([[0.2,0.3,0.5]]) >>> B = ml.matrix([[-1., 0., 0.],[0., -3., 1.],[0., -1., -3.]]) >>> a, A = MonocyclicPHFromME(b, B) >>> print(a) [[ 0.00551 0.00903 0.01694 0.01522 0.00535 0.00874 0.05249 0.22657 0.66016]] >>> print(A) [[-1. 1. 0. 0. 0. 0. 0. 0. 0. ] [ 0. -2.42265 2.42265 0. 0. 0. 0. 0. 0. ] [ 0. 0. -2.42265 2.42265 0. 0. 0. 0. 0. ] [ 0. 0.26232 0. -2.42265 2.16033 0. 0. 0. 0. ] [ 0. 0. 0. 0. -4.2414 4.2414 0. 0. 0. ] [ 0. 0. 0. 0. 0. -4.2414 4.2414 0. 0. ] [ 0. 0. 0. 0. 0. 0. -4.2414 4.2414 0. ] [ 0. 0. 0. 0. 0. 0. 0. -4.2414 4.2414 ] [ 0. 0. 0. 0. 0. 0. 0. 0. -4.2414 ]] >>> b, B = MinimalRepFromME(a, A, "cont") >>> print(b) [[ 0.00551 0.00903 0.01694 0.01522 0.00535 0.00874 0.05249 0.22657 0.66016]] >>> print(B) [[-1. 1. 0. 0. 0. 0. 0. 0. 0. ] [ 0. -2.42265 2.42265 0. 0. 0. 0. 0. 0. ] [ 0. 0. -2.42265 2.42265 0. 0. 0. 0. 0. ] [ 0. 0.26232 0. -2.42265 2.16033 0. 0. 0. 0. ] [ 0. 0. 0. 0. -4.2414 4.2414 0. 0. 0. ] [ 0. 0. 0. 0. 0. -4.2414 4.2414 0. 0. ] [ 0. 0. 0. 0. 0. 0. -4.2414 4.2414 0. ] [ 0. 0. 0. 0. 0. 0. 0. -4.2414 4.2414 ] [ 0. 0. 0. 0. 0. 0. 0. 0. -4.2414 ]] >>> b, B = MinimalRepFromME(a, A, "obs") >>> print(b) [[ 1.00000e+00 6.93889e-18 5.55112e-17]] >>> print(B) [[ -2.83622e+00 3.62221e-02 -2.22045e-16] [ -1.66096e+01 -3.33690e+00 1.60422e+01] [ 1.16435e+00 -5.17236e-02 -8.26876e-01]] >>> Cm = SimilarityMatrix(B, A) >>> err1 = la.norm(B*Cm-Cm*A) >>> err2 = la.norm(b*Cm-a) >>> print(np.max(err1, err2)) 7.41724165831e-15 >>> b, B = MinimalRepFromME(a, A, "obscont") >>> print(b) [[ 1.00000e+00 6.93889e-18 5.55112e-17]] >>> print(B) [[ -2.83622e+00 3.62221e-02 -2.22045e-16] [ -1.66096e+01 -3.33690e+00 1.60422e+01] [ 1.16435e+00 -5.17236e-02 -8.26876e-01]] >>> Cm = SimilarityMatrix(B, A) >>> err1 = la.norm(B*Cm-Cm*A) >>> err2 = la.norm(b*Cm-a) >>> print(np.max(err1, err2)) 7.41724165831e-15 >>> b, B = MinimalRepFromME(a, A, "moment") >>> print(b) [[ 0.33333 0.33333 0.33333]] >>> print(B) [[-2.19048 1.92219 -3.36981] [-1.07693 -2.3906 0.83162] [-0.51037 0.8033 -2.41893]] >>> Cm = SimilarityMatrix(B, A) >>> err1 = la.norm(B*Cm-Cm*A) >>> err2 = la.norm(b*Cm-a) >>> print(np.max(err1, err2)) 2.3964518895e-15