butools.mam.GeneralFluidSolve

butools.mam.GeneralFluidSolve()
Matlab: [mass0, ini, K, clo] = GeneralFluidSolve (Q, R, Q0, prec)
Mathematica: {mass0, ini, K, clo} = GeneralFluidSolve [Q, R, Q0, prec]
Python/Numpy: mass0, ini, K, clo = GeneralFluidSolve (Q, R, Q0, prec)

Returns the parameters of the matrix-exponentially distributed stationary distribution of a general Markovian fluid model, where the fluid rates associated with the states of the background process can be arbitrary (zero is allowed as well).

Using the returned 4 parameters the stationary solution can be obtained as follows.

The probability that the fluid level is zero while being in different states of the background process is given by vector mass0.

The density that the fluid level is x while being in different states of the background process is

\[\pi(x)=ini\cdot e^{K x}\cdot clo.\]
Parameters:

Q : matrix, shape (N,N)

The generator of the background Markov chain

R : diagonal matrix, shape (N,N)

The diagonal matrix of the fluid rates associated with the different states of the background process

Q0 : matrix, shape (N,N), optional

The generator of the background Markov chain at level 0. If not provided, or empty, then Q0=Q is assumed. The default value is empty.

precision : double, optional

Numerical precision for computing the fundamental matrix. The default value is 1e-14
Returns:

mass0 : matrix, shape (1,Np+Nm)

The stationary probability vector of zero level

ini : matrix, shape (1,Np)

The initial vector of the stationary density

K : matrix, shape (Np,Np)

The matrix parameter of the stationary density

clo : matrix, shape (Np,Np+Nm)

The closing matrix of the stationary density

Examples

For Matlab:

>>> Q = [-6., 1., 3., 2., 0., 0.; 6., -10., 2., 0., 2., 0.; 3., 7., -12., 0., 0., 2.; 5., 0., 0., -9., 1., 3.; 0., 5., 0., 6., -13., 2.; 0., 0., 5., 3., 7., -15.];
>>> R = [2., 0., 0., 0., 0., 0.; 0., -4., 0., 0., 0., 0.; 0., 0., -12., 0., 0., 0.; 0., 0., 0., 6., 0., 0.; 0., 0., 0., 0., 0., 0.; 0., 0., 0., 0., 0., -8.];
>>> x = 0.7;
>>> [mass0, ini, K, clo] = GeneralFluidSolve(Q, R);
Final Residual Error for Psi:    7.6328e-16
>>> disp(mass0);
            0     0.082246     0.069492            0     0.023812     0.020724
>>> disp(ini);
      0.70195      0.20505
>>> disp(K);
      -2.4698       1.1349
        1.295      -1.1686
>>> disp(clo);
          0.5     0.061087     0.054574            0      0.01618     0.012595
            0     0.055389     0.043116      0.16667     0.038913     0.032631
>>> pdfAtX = ini*expm(K*x)*clo;
>>> disp(pdfAtX);
      0.12566     0.031849     0.026557     0.049637     0.015655     0.012884

For Mathematica:

>>> Q = {{-6., 1., 3., 2., 0., 0.},{6., -10., 2., 0., 2., 0.},{3., 7., -12., 0., 0., 2.},{5., 0., 0., -9., 1., 3.},{0., 5., 0., 6., -13., 2.},{0., 0., 5., 3., 7., -15.}};
>>> R = {{2., 0., 0., 0., 0., 0.},{0., -4., 0., 0., 0., 0.},{0., 0., -12., 0., 0., 0.},{0., 0., 0., 6., 0., 0.},{0., 0., 0., 0., 0., 0.},{0., 0., 0., 0., 0., -8.}};
>>> x = 0.7;
>>> {mass0, ini, K, clo} = GeneralFluidSolve[Q, R];
"Final Residual Error for Psi: "8.049116928532385*^-16
>>> Print[mass0];
{0., 0.08224612885906893, 0.06949240017690633, 0., 0.023812479297703082, 0.02072428187885746}
>>> Print[ini];
{0.7019539736851326, 0.20504772142279093}
>>> Print[K];
{{-2.4697521505145907, 1.134862595858089},
 {1.2950177373696095, -1.1686305573375777}}
>>> Print[clo];
{{0.5, 0.061087420281185266, 0.05457444259943262, 0., 0.016179786329552276, 0.012594625960258435},
 {0., 0.05538937584186795, 0.04311604966168951, 0.16666666666666666, 0.03891262165047115, 0.03263123758653177}}
>>> pdfAtX = ini.MatrixExponential[K*x].clo;
>>> Print[pdfAtX];
{0.7019539736851326, 0.20504772142279093} . MatrixExponential[{{-1.7288265053602134, 0.7944038171006622}, {0.9065124161587266, -0.8180413901363044}}] . {{0.5, 0.061087420281185266, 0.05457444259943262, 0., 0.016179786329552276, 0.012594625960258435}, {0., 0.05538937584186795, 0.04311604966168951, 0.16666666666666666, 0.03891262165047115, 0.03263123758653177}}

For Python/Numpy:

>>> Q = ml.matrix([[-6., 1., 3., 2., 0., 0.],[6., -10., 2., 0., 2., 0.],[3., 7., -12., 0., 0., 2.],[5., 0., 0., -9., 1., 3.],[0., 5., 0., 6., -13., 2.],[0., 0., 5., 3., 7., -15.]])
>>> R = ml.matrix([[2., 0., 0., 0., 0., 0.],[0., -4., 0., 0., 0., 0.],[0., 0., -12., 0., 0., 0.],[0., 0., 0., 6., 0., 0.],[0., 0., 0., 0., 0., 0.],[0., 0., 0., 0., 0., -8.]])
>>> x = 0.7
>>> mass0, ini, K, clo = GeneralFluidSolve(Q, R)
Final Residual Error for G:  6.661338147750939e-16
>>> print(mass0)
[[ 0.       0.08225  0.06949  0.       0.02381  0.02072]]
>>> print(ini)
[[ 0.70195  0.20505]]
>>> print(K)
[[-2.46975  1.13486]
 [ 1.29502 -1.16863]]
>>> print(clo)
[[ 0.5      0.06109  0.05457  0.       0.01618  0.01259]
 [ 0.       0.05539  0.04312  0.16667  0.03891  0.03263]]
>>> pdfAtX = ini*la.expm(K*x)*clo
>>> print(pdfAtX)
[[ 0.12566  0.03185  0.02656  0.04964  0.01566  0.01288]]