Diff of the two buildlogs:

--
--- b1/build.log	2024-04-28 07:58:04.281442279 +0000
+++ b2/build.log	2024-04-28 08:07:43.824733578 +0000
@@ -1,6 +1,6 @@
 I: pbuilder: network access will be disabled during build
-I: Current time: Sat Apr 27 19:45:12 -12 2024
-I: pbuilder-time-stamp: 1714290312
+I: Current time: Sun Apr 28 21:58:21 +14 2024
+I: pbuilder-time-stamp: 1714291101
 I: Building the build Environment
 I: extracting base tarball [/var/cache/pbuilder/unstable-reproducible-base.tgz]
 I: copying local configuration
@@ -32,52 +32,84 @@
 dpkg-source: info: applying octave7.patch
 I: Not using root during the build.
 I: Installing the build-deps
-I: user script /srv/workspace/pbuilder/7751/tmp/hooks/D02_print_environment starting
+I: user script /srv/workspace/pbuilder/11609/tmp/hooks/D01_modify_environment starting
+debug: Running on virt64z.
+I: Changing host+domainname to test build reproducibility
+I: Adding a custom variable just for the fun of it...
+I: Changing /bin/sh to bash
+'/bin/sh' -> '/bin/bash'
+lrwxrwxrwx 1 root root 9 Apr 28 07:58 /bin/sh -> /bin/bash
+I: Setting pbuilder2's login shell to /bin/bash
+I: Setting pbuilder2's GECOS to second user,second room,second work-phone,second home-phone,second other
+I: user script /srv/workspace/pbuilder/11609/tmp/hooks/D01_modify_environment finished
+I: user script /srv/workspace/pbuilder/11609/tmp/hooks/D02_print_environment starting
 I: set
-  BUILDDIR='/build/reproducible-path'
-  BUILDUSERGECOS='first user,first room,first work-phone,first home-phone,first other'
-  BUILDUSERNAME='pbuilder1'
-  BUILD_ARCH='armhf'
-  DEBIAN_FRONTEND='noninteractive'
-  DEB_BUILD_OPTIONS='buildinfo=+all reproducible=+all parallel=3 '
-  DISTRIBUTION='unstable'
-  HOME='/root'
-  HOST_ARCH='armhf'
+  BASH=/bin/sh
+  BASHOPTS=checkwinsize:cmdhist:complete_fullquote:extquote:force_fignore:globasciiranges:globskipdots:hostcomplete:interactive_comments:patsub_replacement:progcomp:promptvars:sourcepath
+  BASH_ALIASES=()
+  BASH_ARGC=()
+  BASH_ARGV=()
+  BASH_CMDS=()
+  BASH_LINENO=([0]="12" [1]="0")
+  BASH_LOADABLES_PATH=/usr/local/lib/bash:/usr/lib/bash:/opt/local/lib/bash:/usr/pkg/lib/bash:/opt/pkg/lib/bash:.
+  BASH_SOURCE=([0]="/tmp/hooks/D02_print_environment" [1]="/tmp/hooks/D02_print_environment")
+  BASH_VERSINFO=([0]="5" [1]="2" [2]="21" [3]="1" [4]="release" [5]="arm-unknown-linux-gnueabihf")
+  BASH_VERSION='5.2.21(1)-release'
+  BUILDDIR=/build/reproducible-path
+  BUILDUSERGECOS='second user,second room,second work-phone,second home-phone,second other'
+  BUILDUSERNAME=pbuilder2
+  BUILD_ARCH=armhf
+  DEBIAN_FRONTEND=noninteractive
+  DEB_BUILD_OPTIONS='buildinfo=+all reproducible=+all parallel=4 '
+  DIRSTACK=()
+  DISTRIBUTION=unstable
+  EUID=0
+  FUNCNAME=([0]="Echo" [1]="main")
+  GROUPS=()
+  HOME=/root
+  HOSTNAME=i-capture-the-hostname
+  HOSTTYPE=arm
+  HOST_ARCH=armhf
   IFS=' 	
   '
-  INVOCATION_ID='f9100ee514224603be68b3f1a967db0d'
-  LANG='C'
-  LANGUAGE='en_US:en'
-  LC_ALL='C'
-  MAIL='/var/mail/root'
-  OPTIND='1'
-  PATH='/usr/sbin:/usr/bin:/sbin:/bin:/usr/games'
-  PBCURRENTCOMMANDLINEOPERATION='build'
-  PBUILDER_OPERATION='build'
-  PBUILDER_PKGDATADIR='/usr/share/pbuilder'
-  PBUILDER_PKGLIBDIR='/usr/lib/pbuilder'
-  PBUILDER_SYSCONFDIR='/etc'
-  PPID='7751'
-  PS1='# '
-  PS2='> '
+  INVOCATION_ID=0417988d96c64922843e8b396f09ad4f
+  LANG=C
+  LANGUAGE=it_CH:it
+  LC_ALL=C
+  MACHTYPE=arm-unknown-linux-gnueabihf
+  MAIL=/var/mail/root
+  OPTERR=1
+  OPTIND=1
+  OSTYPE=linux-gnueabihf
+  PATH=/usr/sbin:/usr/bin:/sbin:/bin:/usr/games:/i/capture/the/path
+  PBCURRENTCOMMANDLINEOPERATION=build
+  PBUILDER_OPERATION=build
+  PBUILDER_PKGDATADIR=/usr/share/pbuilder
+  PBUILDER_PKGLIBDIR=/usr/lib/pbuilder
+  PBUILDER_SYSCONFDIR=/etc
+  PIPESTATUS=([0]="0")
+  POSIXLY_CORRECT=y
+  PPID=11609
   PS4='+ '
-  PWD='/'
-  SHELL='/bin/bash'
-  SHLVL='2'
-  SUDO_COMMAND='/usr/bin/timeout -k 18.1h 18h /usr/bin/ionice -c 3 /usr/bin/nice /usr/sbin/pbuilder --build --configfile /srv/reproducible-results/rbuild-debian/r-b-build.MjojLgc3/pbuilderrc_6ygR --distribution unstable --hookdir /etc/pbuilder/first-build-hooks --debbuildopts -b --basetgz /var/cache/pbuilder/unstable-reproducible-base.tgz --buildresult /srv/reproducible-results/rbuild-debian/r-b-build.MjojLgc3/b1 --logfile b1/build.log octave-queueing_1.2.7-5.dsc'
-  SUDO_GID='112'
-  SUDO_UID='106'
-  SUDO_USER='jenkins'
-  TERM='unknown'
-  TZ='/usr/share/zoneinfo/Etc/GMT+12'
-  USER='root'
-  _='/usr/bin/systemd-run'
-  http_proxy='http://10.0.0.15:3142/'
+  PWD=/
+  SHELL=/bin/bash
+  SHELLOPTS=braceexpand:errexit:hashall:interactive-comments:posix
+  SHLVL=3
+  SUDO_COMMAND='/usr/bin/timeout -k 24.1h 24h /usr/bin/ionice -c 3 /usr/bin/nice -n 11 /usr/bin/unshare --uts -- /usr/sbin/pbuilder --build --configfile /srv/reproducible-results/rbuild-debian/r-b-build.MjojLgc3/pbuilderrc_N7Tt --distribution unstable --hookdir /etc/pbuilder/rebuild-hooks --debbuildopts -b --basetgz /var/cache/pbuilder/unstable-reproducible-base.tgz --buildresult /srv/reproducible-results/rbuild-debian/r-b-build.MjojLgc3/b2 --logfile b2/build.log octave-queueing_1.2.7-5.dsc'
+  SUDO_GID=110
+  SUDO_UID=103
+  SUDO_USER=jenkins
+  TERM=unknown
+  TZ=/usr/share/zoneinfo/Etc/GMT-14
+  UID=0
+  USER=root
+  _='I: set'
+  http_proxy=http://10.0.0.15:3142/
 I: uname -a
-  Linux virt32b 6.1.0-20-armmp-lpae #1 SMP Debian 6.1.85-1 (2024-04-11) armv7l GNU/Linux
+  Linux i-capture-the-hostname 6.1.0-20-arm64 #1 SMP Debian 6.1.85-1 (2024-04-11) aarch64 GNU/Linux
 I: ls -l /bin
   lrwxrwxrwx 1 root root 7 Apr 25 07:43 /bin -> usr/bin
-I: user script /srv/workspace/pbuilder/7751/tmp/hooks/D02_print_environment finished
+I: user script /srv/workspace/pbuilder/11609/tmp/hooks/D02_print_environment finished
  -> Attempting to satisfy build-dependencies
  -> Creating pbuilder-satisfydepends-dummy package
 Package: pbuilder-satisfydepends-dummy
@@ -678,7 +710,7 @@
 Get: 549 http://deb.debian.org/debian unstable/main armhf texlive all 2023.20240207-1 [18.2 kB]
 Get: 550 http://deb.debian.org/debian unstable/main armhf texlive-font-utils all 2023.20240207-1 [7046 kB]
 Get: 551 http://deb.debian.org/debian unstable/main armhf texlive-plain-generic all 2023.20240207-1 [29.0 MB]
-Fetched 282 MB in 6s (47.9 MB/s)
+Fetched 282 MB in 13s (21.0 MB/s)
 debconf: delaying package configuration, since apt-utils is not installed
 Selecting previously unselected package readline-common.
 (Reading database ... 
(Reading database ... 5%
(Reading database ... 10%
(Reading database ... 15%
(Reading database ... 20%
(Reading database ... 25%
(Reading database ... 30%
(Reading database ... 35%
(Reading database ... 40%
(Reading database ... 45%
(Reading database ... 50%
(Reading database ... 55%
(Reading database ... 60%
(Reading database ... 65%
(Reading database ... 70%
(Reading database ... 75%
(Reading database ... 80%
(Reading database ... 85%
(Reading database ... 90%
(Reading database ... 95%
(Reading database ... 100%
(Reading database ... 19443 files and directories currently installed.)
@@ -2954,7 +2986,11 @@
 Building tag database...
  -> Finished parsing the build-deps
 I: Building the package
-I: Running cd /build/reproducible-path/octave-queueing-1.2.7/ && env PATH="/usr/sbin:/usr/bin:/sbin:/bin:/usr/games" HOME="/nonexistent/first-build" dpkg-buildpackage -us -uc -b && env PATH="/usr/sbin:/usr/bin:/sbin:/bin:/usr/games" HOME="/nonexistent/first-build" dpkg-genchanges -S  > ../octave-queueing_1.2.7-5_source.changes
+I: user script /srv/workspace/pbuilder/11609/tmp/hooks/A99_set_merged_usr starting
+Not re-configuring usrmerge for unstable
+I: user script /srv/workspace/pbuilder/11609/tmp/hooks/A99_set_merged_usr finished
+hostname: Name or service not known
+I: Running cd /build/reproducible-path/octave-queueing-1.2.7/ && env PATH="/usr/sbin:/usr/bin:/sbin:/bin:/usr/games:/i/capture/the/path" HOME="/nonexistent/second-build" dpkg-buildpackage -us -uc -b && env PATH="/usr/sbin:/usr/bin:/sbin:/bin:/usr/games:/i/capture/the/path" HOME="/nonexistent/second-build" dpkg-genchanges -S  > ../octave-queueing_1.2.7-5_source.changes
 dpkg-buildpackage: info: source package octave-queueing
 dpkg-buildpackage: info: source version 1.2.7-5
 dpkg-buildpackage: info: source distribution unstable
@@ -3642,67 +3678,195 @@
 Checking package...
 Run the unit tests...
 Checking m files ...
-[inst/ctmcfpt.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/ctmcfpt.m
+[inst/qncsmvaap.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsmvaap.m
+***** test
+ fail( "qncsmvaap()", "Invalid" );
+ fail( "qncsmvaap( 10, [1 2], [1 2 3] )", "S, V and m" );
+ fail( "qncsmvaap( 10, [-1 1], [1 1] )", ">= 0" );
+ fail( "qncsmvaap( 10, [1 2], [1 2], [1 2] )", "supports");
+ fail( "qncsmvaap( 10, [1 2], [1 2], [1 1], 0, -1)", "tol");
+***** test
+ # Example p. 117 Lazowska et al.
+ S = [0.605 2.1 1.35];
+ V = [1 1 1];
+ N = 3;
+ Z = 15;
+ m = 1;
+ [U R Q X] = qncsmvaap(N, S, V, m, Z);
+ Rs = dot(V,R);
+ Xs = N/(Z+Rs);
+ assert( Q, [0.0973 0.4021 0.2359], 1e-3 );
+ assert( Xs, 0.1510, 1e-3 );
+ assert( Rs, 4.87, 1e-3 );
 ***** demo
- Q = [ -1.0  0.9  0.1; ...
-        0.1 -1.0  0.9; ...
-        0.9  0.1 -1.0 ];
- M = ctmcfpt(Q)
- m = ctmcfpt(Q,1,3)
+ S = [ 0.125 0.3 0.2 ];
+ V = [ 16 10 5 ];
+ N = 30;
+ m = ones(1,3);
+ Z = 4;
+ Xmva = Xapp = Rmva = Rapp = zeros(1,N);
+ for n=1:N
+   [U R Q X] = qncsmva(n,S,V,m,Z);
+   Xmva(n) = X(1)/V(1);
+   Rmva(n) = dot(R,V);
+   [U R Q X] = qncsmvaap(n,S,V,m,Z);
+   Xapp(n) = X(1)/V(1);
+   Rapp(n) = dot(R,V);
+ endfor
+ subplot(2,1,1);
+ plot(1:N, Xmva, ";Exact;", "linewidth", 2, 1:N, Xapp, "x;Approximate;", "markersize", 7);
+ legend("location","southeast"); legend("boxoff");
+ ylabel("Throughput X(n)");
+ subplot(2,1,2);
+ plot(1:N, Rmva, ";Exact;", "linewidth", 2, 1:N, Rapp, "x;Approximate;", "markersize", 7);
+ legend("location","southeast"); legend("boxoff");
+ ylabel("Response Time R(n)");
+ xlabel("Number of Requests n");
+2 tests, 2 passed, 0 known failure, 0 skipped
+[inst/ctmcbd.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/ctmcbd.m
 ***** test
- N = 10;
- Q = reshape(1:N^2,N,N);
- Q(1:N+1:end) = 0;
- Q -= diag(sum(Q,2));
- M = ctmcfpt(Q);
- assert( all(diag(M) < 10*eps) );
+ birth = [ 1 1 1 ];
+ death = [ 2 2 2 ];
+ Q = ctmcbd( birth, death );
+ assert( ctmc(Q), [ 8/15 4/15 2/15 1/15 ], 1e-5 );
 1 test, 1 passed, 0 known failure, 0 skipped
-[inst/qnosvisits.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnosvisits.m
+[inst/dtmctaexps.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/dtmctaexps.m
 ***** test
- fail( "qnosvisits([0 .5; .5 0],[0 -1])", "contains negative" );
- fail( "qnosvisits([1 1 1; 1 1 1], [1 1])", "square" );
+ P = dtmcbd([1 1 1 1], [0 0 0 0]);
+ p0 = [1 0 0 0 0];
+ L = dtmctaexps(P,p0);
+ assert( L, [.25 .25 .25 .25 0], 10*eps );
+1 test, 1 passed, 0 known failure, 0 skipped
+[inst/qnsolve.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnsolve.m
 ***** test
+ # Example 8.7 p. 349 Bolch et al.
+ N = 3;
+ Q1 = qnmknode( "m/m/m-fcfs", .5, 2 );
+ Q2 = qnmknode( "m/m/m-fcfs", 1/1.667 );
+ Q3 = qnmknode( "m/m/m-fcfs", 1/1.25 );
+ Q4 = qnmknode( "m/m/m-fcfs", 1./[1 2 3] );
+ V = [ 1 .5 .5 1 ];
+ [U R Q X] = qnsolve("closed",N, { Q1, Q2, Q3, Q4 }, V);
+ assert( Q, [0.624 0.473 0.686 1.217], 1e-3 );
+ assert( R, [0.512 0.776 1.127 1], 1e-3 );
+warning: the '.+' operator was deprecated in version 7 and will not be allowed in a future version of Octave; please use '+' instead; near line 424 of file '/build/reproducible-path/octave-queueing-1.2.7/debian/octave-queueing/usr/share/octave/packages/queueing-1.2.7/qnsolve.m'
+warning: called from
+    __test__ at line 10 column 15
+    test at line 682 column 11
+    /tmp/tmp.qG5gOSRbSf at line 38 column 31
 
- ## Open, single class network
-
- P = [0 0.2 0.5; 1 0 0; 1 0 0];
- lambda = [ 0.1 0.3 0.2 ];
- V = qnosvisits(P,lambda);
- assert( V*P+lambda/sum(lambda),V,1e-5 );
-***** demo
- p = 0.3;
- lambda = 1.2
- P = [0 0.3 0.5; ...
-      1 0   0  ; ...
-      1 0   0  ];
- V = qnosvisits(P,[1.2 0 0])
-***** demo
- P = [ 0 0.4 0.6 0; ...
-       0.2 0 0.2 0.6; ...
-       0 0 0 1; ...
-       0 0 0 0 ];
- lambda = [0.1 0 0 0.3];
- V = qnosvisits(P,lambda);
- S = [2 1 2 1.8];
- m = [3 1 1 2];
- [U R Q X] = qnos( sum(lambda), S, V, m )
-2 tests, 2 passed, 0 known failure, 0 skipped
-[inst/qnosaba.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnosaba.m
 ***** test
- fail( "qnosaba( 0.1, [] )", "vector" );
- fail( "qnosaba( 0.1, [0 -1])", "nonnegative" );
- fail( "qnosaba( 0, [1 2] )", "lambda" );
- fail( "qnosaba( -1, [1 2])", "lambda" );
- fail( "qnosaba( 1, [1 2 3], [1 2] )", "incompatible size");
- fail( "qnosaba( 1, [1 2 3], [-1 2 3] )", "nonnegative");
+ # Example 8.7 p. 349 Bolch et al.
+ N = 3;
+ Q1 = qnmknode( "m/m/m-fcfs", 1/2, 2 );
+ Q2 = qnmknode( "m/m/m-fcfs", 1/1.667 );
+ Q3 = qnmknode( "m/m/m-fcfs", 1/1.25 );
+ Q4 = qnmknode( "-/g/inf", 1 );
+ V = [ 1 .5 .5 1 ];
+ [U R Q X] = qnsolve("closed",N, { Q1, Q2, Q3, Q4 }, V);
+ assert( Q, [0.624 0.473 0.686 1.217], 1e-3 );
+ assert( R, [0.512 0.776 1.127 1], 1e-3 );
 ***** test
- [Xl Xu Rl Ru] = qnosaba( 1, [1 1] );
- assert( Xl, 0 );
- assert( Ru, +inf );
-2 tests, 2 passed, 0 known failure, 0 skipped
+ # Example 8.4 p. 333 Bolch et al.
+ N = 3;
+ Q1 = qnmknode( "m/m/m-fcfs", .5, 2 );
+ Q2 = qnmknode( "m/m/m-fcfs", .6 );
+ Q3 = qnmknode( "m/m/m-fcfs", .8 );
+ Q4 = qnmknode( "-/g/inf", 1 );
+ V = [ 1 .5 .5 1 ];
+ [U R Q X] = qnsolve("closed",N, { Q1, Q2, Q3, Q4 }, V);
+ assert( U(1:3), [.304 .365 .487], 1e-3 );
+ assert( X, [1.218 0.609 0.609 1.218], 1e-3 );
+***** test
+ # Same as above, with center 1 replaced with a load-dependent service center
+ N = 3;
+ Q1 = qnmknode( "m/m/m-fcfs", [.5 .25 .25] );
+ Q2 = qnmknode( "m/m/m-fcfs", .6 );
+ Q3 = qnmknode( "m/m/m-fcfs", .8 );
+ Q4 = qnmknode( "m/m/m-fcfs", [1 1/2 1/3] );
+ V = [ 1 .5 .5 1 ];
+ [U R Q X] = qnsolve("closed",N, { Q1, Q2, Q3, Q4 }, V);
+ assert( U(2:3), [.365 .487], 1e-3 ); ## NOTE: Utilization U(1) is computed differently from M/M/m nodes and load-dependent M/M/1 nodes
+ assert( X, [1.218 0.609 0.609 1.218], 1e-3 );
+***** test
+ # Example 7.4 p. 287 Bolch et al.
+ QQ = { qnmknode( "m/m/m-fcfs", 0.04 ), ...
+        qnmknode( "m/m/m-fcfs", 0.03 ), ...
+        qnmknode( "m/m/m-fcfs", 0.06 ), ...
+        qnmknode( "m/m/m-fcfs", 0.05 ) };
+ P = [ 0 0.5 0.5 0; 1 0 0 0; 0.6 0 0 0; 1 0 0 0 ];
+ lambda = [0 0 0 4];
+ [U R Q X] = qnsolve("open", sum(lambda), QQ, qnosvisits(P,lambda) );
+ assert( X, [20 10 10 4], 1e-4 );
+ assert( U, [0.8 0.3 0.6 0.2], 1e-2 );
+ assert( R, [0.2 0.043 0.15 0.0625], 1e-3 );
+ assert( Q, [4, 0.429 1.5 0.25], 1e-3 );
+***** test
+ V = [1 1; 1 1];
+ Q1 = qnmknode( "m/m/m-fcfs", [1;2] );
+ Q2 = qnmknode( "m/m/m-fcfs", [3;4] );
+ lambda = [3/19 2/19];
+ [U R Q] = qnsolve("open", lambda, { Q1, Q2 }, diag( lambda / sum(lambda) ) * V);
+ assert( U(1,1), 3/19, 1e-6 );
+ assert( U(2,1), 4/19, 1e-6 );
+ assert( R(1,1), 19/12, 1e-6 );
+ assert( R(1,2), 57/2, 1e-6 );
+ assert( Q(1,1), .25, 1e-6 );
+***** test
+ QQ = { qnmknode( "m/m/m-fcfs", [0.2; 0.2], 2 ), ...
+        qnmknode( "-/g/1-ps", [0.4; 0.6] ), ...
+        qnmknode( "-/g/inf", [1; 2] ) };
+ V = [ 1 0.6 0.4; 1 0.3 0.7 ];
+ N = [ 2 1 ];
+ [U R Q X] = qnsolve( "closed", N, QQ, V );
+ assert( Q, [ 0.428 0.726 0.845; 0.108 0.158 0.734 ], 1e-3 );
+ assert( X(1,1), 2.113, 1e-3 ); # CHECK
+ assert( X(2,1), 0.524, 1e-3 ); # CHECK
+ assert( all( all(U(:,[1,2])<=1) ) );
+***** test
+ QQ = { qnmknode( "m/m/m-fcfs", [0.2 0.1 0.1; 0.2 0.1 0.1] ), ...
+        qnmknode( "-/g/1-ps", [0.4; 0.6] ), ...
+        qnmknode( "-/g/inf", [1; 2] ) };
+ V = [ 1 0.6 0.4; 1 0.3 0.7 ];
+ N = [ 2 1 ];
+ [U R Q X] = qnsolve( "closed", N, QQ, V );
+ assert( Q, [ 0.428 0.726 0.845; 0.108 0.158 0.734 ], 1e-3 );
+ assert( X(1,1), 2.113, 1e-3 ); # CHECK
+ assert( X(2,1), 0.524, 1e-3 ); # CHECK
+ assert( all( all(U(:,[1,2])<=1) ) );
+***** test
+ # example p. 26 Schwetman
+ QQ = { qnmknode( "m/m/m-fcfs", [.25; .25] ),
+        qnmknode( "-/g/1-ps", [0; .1] ) };
+ V = [1 0; 1 1];
+ lambda = [1 0];
+ N = [0 3];
+ [U R Q X] = qnsolve( "mixed", lambda, N, QQ, V );
+ assert( U(1,1), .25, 1e-3 );
+ assert( X(1,1), 1.0, 1e-3 );
+ assert( [R(1,1) R(2,1) R(2,2)], [1.201 0.885 0.135], 1e-3 );
+warning: the '.+' operator was deprecated in version 7 and will not be allowed in a future version of Octave; please use '+' instead; near line 430 of file '/build/reproducible-path/octave-queueing-1.2.7/debian/octave-queueing/usr/share/octave/packages/queueing-1.2.7/qncmmva.m'
+warning: called from
+    qnmix at line 187 column 24
+    qnsolve>__qnsolve_mixed at line 514 column 16
+    qnsolve at line 121 column 18
+    __test__ at line 9 column 15
+    test at line 682 column 11
+    /tmp/tmp.qG5gOSRbSf at line 38 column 31
+
+***** demo
+ QQ = { qnmknode( "m/m/m-fcfs", [0.2 0.1 0.1; 0.2 0.1 0.1] ), ...
+        qnmknode( "-/g/1-ps", [0.4; 0.6] ), ...
+        qnmknode( "-/g/inf", [1; 2] ) };
+ V = [ 1 0.6 0.4; ...
+       1 0.3 0.7 ];
+ N = [ 2 1 ];
+ [U R Q X] = qnsolve( "closed", N, QQ, V );
+9 tests, 9 passed, 0 known failure, 0 skipped
 [inst/dtmcmtta.m]
 >>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/dtmcmtta.m
 ***** test
@@ -3803,71 +3967,108 @@
  t = dtmcmtta(P);
  assert( t, [6 5 6 5 0 5 6 5 6], 10*eps );
 5 tests, 5 passed, 0 known failure, 0 skipped
-[inst/qncsmvaap.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsmvaap.m
+[inst/qnmarkov.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnmarkov.m
 ***** test
- fail( "qncsmvaap()", "Invalid" );
- fail( "qncsmvaap( 10, [1 2], [1 2 3] )", "S, V and m" );
- fail( "qncsmvaap( 10, [-1 1], [1 1] )", ">= 0" );
- fail( "qncsmvaap( 10, [1 2], [1 2], [1 2] )", "supports");
- fail( "qncsmvaap( 10, [1 2], [1 2], [1 1], 0, -1)", "tol");
+ S = [5 2.5];
+ P = [0 1; 1 0];
+ C = [3 3];
+ m = [1 1];
+ [U R Q X] = qnmarkov( 3, S, C, P, m );
+ assert( U, [0.9333 0.4667], 1e-4 );
+ assert( X, [0.1867 0.1867], 1e-4 );
+ assert( R, [12.1429 3.9286], 1e-4 );
 ***** test
- # Example p. 117 Lazowska et al.
- S = [0.605 2.1 1.35];
- V = [1 1 1];
- N = 3;
- Z = 15;
- m = 1;
- [U R Q X] = qncsmvaap(N, S, V, m, Z);
- Rs = dot(V,R);
- Xs = N/(Z+Rs);
- assert( Q, [0.0973 0.4021 0.2359], 1e-3 );
- assert( Xs, 0.1510, 1e-3 );
- assert( Rs, 4.87, 1e-3 );
-***** demo
- S = [ 0.125 0.3 0.2 ];
- V = [ 16 10 5 ];
- N = 30;
- m = ones(1,3);
- Z = 4;
- Xmva = Xapp = Rmva = Rapp = zeros(1,N);
- for n=1:N
-   [U R Q X] = qncsmva(n,S,V,m,Z);
-   Xmva(n) = X(1)/V(1);
-   Rmva(n) = dot(R,V);
-   [U R Q X] = qncsmvaap(n,S,V,m,Z);
-   Xapp(n) = X(1)/V(1);
-   Rapp(n) = dot(R,V);
- endfor
- subplot(2,1,1);
- plot(1:N, Xmva, ";Exact;", "linewidth", 2, 1:N, Xapp, "x;Approximate;", "markersize", 7);
- legend("location","southeast"); legend("boxoff");
- ylabel("Throughput X(n)");
- subplot(2,1,2);
- plot(1:N, Rmva, ";Exact;", "linewidth", 2, 1:N, Rapp, "x;Approximate;", "markersize", 7);
- legend("location","southeast"); legend("boxoff");
- ylabel("Response Time R(n)");
- xlabel("Number of Requests n");
-2 tests, 2 passed, 0 known failure, 0 skipped
-[inst/qsmminf.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qsmminf.m
+ S = [1/0.8 1/0.6 1/0.4];
+ P = [0.6 0.3 0.1; 0.2 0.3 0.5; 0.4 0.1 0.5];
+ C = [3 3 3];
+ [U R Q X] = qnmarkov( 3, S, C, P );
+ assert( U, [0.543 0.386 0.797], 1e-3 );
+ assert( Q, [0.873 0.541 1.585], 1e-3 );
+***** xtest
+ S = [2 0.9];
+ C = [7 5];
+ P = [0 1; 1 0];
+ [U R Q X] = qnmarkov( 10, S, C, P );
+ assert( Q, [6.73 3.27], 1e-3 );
+!!!!! known failure
+ASSERT errors for:  assert (Q,[6.73, 3.27],1e-3)
+
+  Location  |  Observed  |  Expected  |  Reason
+    (1)         6.4826        6.73       Abs err 0.2474 exceeds tol 0.001 by 0.2
+    (2)         3.5174        3.27       Abs err 0.2474 exceeds tol 0.001 by 0.2
 ***** test
- fail( "qsmminf( [1 2], [1 2 3] )", "incompatible size");
- fail( "qsmminf( [-1 -1], [1 1] )", ">0" );
-***** demo
- ## Given a M/M/inf and M/M/m queue, compute the steady-state probability pk
- ## of having k requests in the systen.
- lambda = 5;
- mu = 1.1;
- m = 5;
- k = 0:20;
- pk_inf = qsmminf(lambda, mu, k);
- pk_m = qsmmm(lambda, mu, 5, k);
- plot(k, pk_inf, "-o;M/M/\\infty;", "linewidth", 2, ...
-      k, pk_m, "-x;M/M/5;", "linewidth", 2);
- xlabel("N. of requests (k)");
- ylabel("P_k");
- title(sprintf("M/M/\\infty and M/M/%d systems, \\lambda = %g, \\mu = %g", m, lambda, mu));
+ S = [1/0.8 1/0.6 1/0.4];
+ P = [(1-0.667-0.2) 0.667 0.2; 1 0 0; 1 0 0];
+ m = [2 3 1];
+ C = [3 3 3];
+ [U R Q X] = qnmarkov( 3, S, C, P, m );
+ assert( U, [0.590 0.350 0.473], 1e-3 );
+ assert( Q(1:2), [1.290 1.050], 1e-3 );
+***** test
+ p = 0.5; # transition prob. from S2 to S1
+ mu = [1 2]; # Service rates
+ C = [2 1]; # Capacities
+ lambda = [0.5 0]; # arrival rate at service center 1
+
+ PP = [ 0 1; p 0 ];
+ [U R Q X] = qnmarkov( lambda, 1./mu, C, PP );
+ ## Now we generate explicitly the infinitesimal generator matrix
+ ## of the underlying MC.
+ ##    00  01  10  11  20  21
+ QQ = [ 0 0 lambda(1) 0 0 0; ... ## 00
+        mu(2)*(1-p) 0 mu(2)*p lambda(1) 0 0; ... ## 01
+        0 mu(1) 0 0 lambda(1) 0; ... ## 10
+        0 0 mu(2)*(1-p) 0 mu(2)*p lambda(1); ... ## 11
+        0 0 0 mu(1) 0 0; ... ## 20
+        0 0 0 0 mu(2)*(1-p) 0 ]; ## 21
+ ## Complete matrix
+ sum_el = sum(QQ,2);
+ QQ -= diag(sum_el);
+ q = ctmc(QQ);
+ ## Compare results
+ assert( U(1), 1-sum(q([1, 2])), 1e-5 );
+ assert( U(2), 1-sum(q([1,3,5])), 1e-5 );
+***** test
+ P = [0 0.5 0.5; 1 0 0; 1 0 0];
+ C = [6 6 6];
+ S = [1 0.8 1.8];
+ N = 6;
+ [U1 R1 Q1 X1] = qnclosed( N, S, qncsvisits(P) );
+ [U2 R2 Q2 X2] = qnmarkov( N, S, C, P );
+ assert( U1, U2, 1e-6 );
+ assert( R1, R2, 1e-6 );
+ assert( Q1, Q2, 1e-6 );
+ assert( X1, X2, 1e-6 );
+6 tests, 5 passed, 1 known failure, 0 skipped
+[inst/engset.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/engset.m
+***** test
+ fail("erlangb(1, -1)", "positive");
+ fail("erlangb(-1, 1)", "positive");
+ fail("erlangb(1, 0)", "positive");
+ fail("erlangb(0, 1)", "positive");
+ fail("erlangb('foo',1)", "positive");
+ fail("erlangb(1,'bar')", "positive");
+ fail("erlangb([1 1],[1 1 1])","common size");
+1 test, 1 passed, 0 known failure, 0 skipped
+[inst/qncsbsb.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsbsb.m
+***** test
+ fail("qncsbsb(-1,0)", "N must be");
+ fail("qncsbsb(1,[])", "nonempty");
+ fail("qncsbsb(1,[-1 2])", "nonnegative");
+ fail("qncsbsb(1,[1 2],[1 2 3])", "incompatible size");
+ fail("qncsbsb(1,[1 2 3],[1 2 3],[1 2])", "incompatible size");
+ fail("qncsbsb(1,[1 2 3],[1 2 3],[1 2 1])", "M/M/1 servers");
+ fail("qncsbsb(1,[1 2 3],[1 2 3],[1 1 1],-1)", "nonnegative");
+ fail("qncsbsb(1,[1 2 3],[1 2 3],[1 1 1],[0 0])", "scalar");
+1 test, 1 passed, 0 known failure, 0 skipped
+[inst/erlangc.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/erlangc.m
+***** test
+ fail("erlangc('foo',1)", "positive");
+ fail("erlangc(1,'bar')", "positive");
 1 test, 1 passed, 0 known failure, 0 skipped
 [inst/qncmbsb.m]
 >>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncmbsb.m
@@ -3887,7 +4088,7 @@
     fail at line 130 column 7
     __test__ at line 3 column 2
     test at line 682 column 11
-    /tmp/tmp.cnygCL2JrV at line 62 column 31
+    /tmp/tmp.qG5gOSRbSf at line 86 column 31
 
 ***** test
  [Xl Xu Rl Ru] = qncmbsb([0 0], [1 2 3; 1 2 3]);
@@ -3924,643 +4125,6 @@
  ylim([0, 700]);
  title("Class 2 response time"); legend("location", "northwest"); legend("boxoff");
 2 tests, 2 passed, 0 known failure, 0 skipped
-[inst/qncmmva.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncmmva.m
-***** test
- S = [1 1 2; 1 1 1];
- V = [1 1 1; 1 1 1];
- N = [1 1];
- m = [1 1 2];
- fail( "qncmmva(N)" );
- fail( "qncmmva(N,S,V,m)", "independent" );
- S = [0 0 0; 1 1 1];
- fail( "qncmmva(N,S,V,m)", "must contain at least" );
- S = [1 2 3; 1 2 3];
- N = [1 1];
- V = zeros(3,2,3);
- fail( "qncmmva(N,S,V)", "size mismatch" );
- fail( "qncmmva([0.3 1], [1 2; 3 4])", "integer");
- fail( "qncmmva([-1 0], [1 2; 3 4])", "nonnegative");
-warning: the '.+' operator was deprecated in version 7 and will not be allowed in a future version of Octave; please use '+' instead; near line 430 of file '/build/reproducible-path/octave-queueing-1.2.7/debian/octave-queueing/usr/share/octave/packages/queueing-1.2.7/qncmmva.m'
-warning: called from
-    fail at line 130 column 7
-    __test__ at line 7 column 2
-    test at line 682 column 11
-    /tmp/tmp.cnygCL2JrV at line 70 column 31
-
-warning: the '.+' operator was deprecated in version 7 and will not be allowed in a future version of Octave; please use '+' instead; near line 430 of file '/build/reproducible-path/octave-queueing-1.2.7/debian/octave-queueing/usr/share/octave/packages/queueing-1.2.7/qncmmva.m'
-warning: called from
-    print_usage at line 62 column 18
-    qncmmva at line 204 column 5
-    fail at line 130 column 7
-    __test__ at line 7 column 2
-    test at line 682 column 11
-    /tmp/tmp.cnygCL2JrV at line 70 column 31
-
-***** test
- S = [1 1 1; 1 1 1];
- N = [0 0];
- [U R Q X] = qncmmva(N, S);
- assert( U, 0*S );
- assert( R, 0*S );
- assert( Q, 0*S );
- assert( X, 0*S );
-***** test
- S = [1 1 1; 1 1 1];
- V = [1 1 1; 1 1 1];
- N = [0 0];
- m = [2 2 2];
- [U R Q X] = qncmmva(N, S, V, m);
- assert( U, 0*S );
- assert( R, 0*S );
- assert( Q, 0*S );
- assert( X, 0*S );
-***** test
- S = [ 1/10 1/3; 2/5 1 ];
- V = [ 10 9; 5 4 ];
- N = [ 1 1 ];
- [U R Q X] = qncmmva(N,S,V);
- assert( Q, [ 4/19 15/19; 5/19 14/19 ], 1e-3 );
- assert( R .* V, [ 4/3 5; 5/2 7 ], 1e-3 );
- assert( diag( X ./ V )', [ 3/19 2/19 ], 1e-3 );
- assert( all(U(:)<=1) );
- assert( Q, R.*X, 1e-5 ); # Little's Law
-***** test
- S = [0.02 0.2 0.4 0.6];
- V = [1 0.4 0.2 0.1];
- N = [6];
- [U R Q X] = qncmmva( N, S, V );
- assert( Q, [0.244 2.261 2.261 1.234], 1e-3 );
- assert( R, [0.025 0.570 1.140 1.244], 1e-3 );
- assert( X, [9.920 3.968 1.984 0.992], 1e-3 );
- assert( U, [0.198 0.794 0.794 0.595], 1e-3 );
-***** test
- S = [1 0 .025; 0 15 .5];
- V = [1 0 1; 0 1 1];
- N = [2 1];
- m = [-1 -1 1];
- [U R Q X] = qncmmva(N,S,V,m);
- assert( R(1,1), 1, 1e-3 );
- assert( R(2,2), 15, 1e-3 );
- assert( R(1,3), .027, 1e-3 );
- assert( R(2,3), .525, 1e-3  );
- assert( X(1,1)+X(1,2), 1.949, 1e-3 );
- assert( X(2,1)+X(2,2), 0.064, 1e-3 );
- assert( sum(Q,1), [1.949, .966, .085], 1e-3 );
- assert( all(U(:,3)<=1) );
- assert( Q, R.*X, 1e-5 ); # Little's Law
-***** test
- S = [1 0 .025; 0 15 .5];
- V = [1 0 1; 0 1 1];
- N = [15 5];
- m = [-1 -1 1];
- [U R Q X] = qncmmva(N,S,V,m);
- # FIXME: I replaced 14.3->14.323
- assert( U, [14.323 0 .358; 0 4.707 .157], 1e-3 );
- # FIXME: I replaced 14.3->14.323
- assert( X, [14.323 0 14.323; 0 .314 .314 ], 1e-3 );
- # FIXME: I replaced 14.3->14.323
- assert( Q, [14.323 0 .677; 0 4.707 .293 ], 1e-3 );
- assert( R, [1 0 .047; 0 15 .934 ], 1e-3 );
- assert( Q, R.*X, 1e-5 ); # Little's Law
-***** test
- S = [ 0.2 0.4 1; 0.2 0.6 2 ];
- V = [ 1 0.6 0.4; 1 0.3 0.7 ];
- N = [ 2 1 ];
- m = [ 2 1 -1 ];
- [U R Q X] = qncmmva(N,S,V,m);
- assert( Q, [ 0.428 0.726 0.845; 0.108 0.158 0.734 ], 1e-3 );
- assert( X(1,1), 2.113, 1e-3 ); # CHECK
- assert( X(2,1), 0.524, 1e-3 ); # CHECK
- assert( all(U(:)<=1) );
- assert( Q, R.*X, 1e-5 ); # Little's Law
-***** test
- C = 2; # two classes
- K = 4; # four servers
- S = V = zeros(C,K);
- S(1,:) = linspace(1,2,K);
- S(2,:) = linspace(2,3,K);
- V(1,:) = linspace(4,1,K);
- V(2,:) = linspace(6,3,K);
- N = [10 0]; # class 2 has no customers
- [U1 R1 Q1 X1] = qncmmva(N,S,V);
- [U2 R2 Q2 X2] = qncsmva(N(1),S(1,:),V(1,:));
- assert( U1(1,:), U2, 1e-5 );
- assert( R1(1,:), R2, 1e-5 );
- assert( Q1(1,:), Q2, 1e-5 );
- assert( X1(1,:), X2, 1e-5 );
-***** test
- Z = [1 15];
- V = [1; 1];
- S = [.025; .5];
- N = [15; 5];
- [U R Q X] = qncmmva(N, S, V, 1, Z);
- assert( U, [.358; .157], 1e-3 );
- assert( Q, [.677; .293], 1e-3 );
- assert( X, [14.323; .314], 1e-3 ); ## NOTE: X(1,1) = 14.3 in Schwetman
- assert( R, [.047; .934], 1e-3 );
- assert( Q, R.*X, 1e-5 ); # Little's Law
-***** test
- C = 2;
- K = 6;
- P = zeros(C,K,C,K);
- P(1,1,1,2) = P(2,1,2,2) = 1;
- P(1,2,1,3) = P(1,2,1,4) = P(1,2,1,5) = P(1,2,1,6) = .25;
- P(2,2,2,3) = P(2,2,2,4) = P(2,2,2,5) = P(2,2,2,6) = .25;
- P(1,3,1,1) = P(1,4,1,1) = P(1,5,1,1) = P(1,6,1,1) = .9;
- P(1,3,1,2) = P(1,4,1,2) = P(1,5,1,2) = P(1,6,1,2) = .1;
- P(2,3,2,1) = P(2,4,2,1) = P(2,5,2,1) = P(2,6,2,1) = .05;
- P(2,3,2,2) = P(2,4,2,2) = P(2,5,2,2) = P(2,6,2,2) = .95;
- N = [40 4];
- S = [ 5.0 .010 .035 .035 .035 .035; ...
-      10.0 .100 .035 .035 .035 .035 ];
- V = qncmvisits(P);
- [U R Q X] = qncmmva(N, S, V, [-1 1 1 1 1 1]);
- # FIXME: The results below were computed with JMVA; the numbers
- # in the paper are different (wrong?!?)!!
- assert( U, [39.457941 0.087684 0.076724 0.076724 0.076724 0.076724; ...
-              2.772704 0.554541 0.048522 0.048522 0.048522 0.048522 ], 1e-5 );
- assert( R.*V, [5 0.024363 0.011081 0.011081 0.011081 0.011081; ...
-                10 3.636155 0.197549 0.197549 0.197549 0.197549 ], 1e-5 );
- assert( Q(:,1), [39.457941 2.772704]', 1e-5 );
- assert( Q(:,2), [0.192262 1.008198]', 1e-5 );
- assert( Q(:,3), [0.087449 0.054775]', 1e-5 );
- assert( Q(:,4), Q(:,5), 1e-5 );
- assert( Q(:,5), Q(:,6), 1e-5 );
- assert( X(:,1), [7.891588 0.277270]', 1e-5 );
- assert( X(:,2), [8.768431 5.545407]', 1e-5 );
- assert( X(:,3), [2.192108 1.386352]', 1e-5 );
- assert( X(:,4), X(:,5), 1e-5 );
- assert( X(:,5), X(:,6), 1e-5 );
- assert( Q, R.*X, 1e-5 ); # Little's Law
-***** test
- C = 2; # two classes
- K = 4; # four servers
- C = 2; K = 4;
- P = zeros(C,K,C,K);
- S = zeros(C,K);
-
- # Routing
-
- # class 1 routing
- P(1,1,1,1) = .05;
- P(1,1,1,2) = .45;
- P(1,1,1,3) = .5;
- P(1,2,1,1) = 1;
- P(1,3,1,1) = 1;
- # class 2 routing
- P(2,1,2,1) = .01;
- P(2,1,2,3) = .5;
- P(2,1,2,4) = .49;
- P(2,3,2,1) = 1;
- P(2,4,2,1) = 1;
-
- # Compute visits
-
- V = qncmvisits(P);
-
- # Define population and service times
- 
- N = [3 2];
- S = [0.01 0.09 0.10 0.08; ...
-      0.05 0.09 0.10 0.08];
- [U1 R1 Q1 X1] = qncmmva(N,S,V); # this invokes __qncmmva_nocs
- [U2 R2 Q2 X2] = qncmmva(N,S,P); # this invokes __qncmmva_cs
- assert( U2, U1, 1e-5 );
- assert( R2, R1, 1e-5 );
- assert( Q2, Q1, 1e-5 );
- assert( X2, X1, 1e-5 );
-***** test
- S = [1 0 .025; 0 15 .5];
- V = [1 0 1; 0 1 1];
- N = [15 5];
- m = [-1 -1 1];
- [U R Q X] = qncmmva(N,S,V,m);
- assert( U, [14.323 0 .358; 0 4.707 .157], 1e-3 );
- assert( R, [1.0 0 .047; 0 15 .934], 1e-3 );
- assert( Q, [14.323 0 .677; 0 4.707 .293], 1e-3 );
- assert( X, [14.323 0 14.323; 0 .314 .314], 1e-3 );
-***** test
- S = [.025 1 15; .5 1 15 ];
- P = zeros(2,3,2,3);
- P(1,1,1,2) = P(1,2,1,1) = 1;
- P(2,1,2,3) = P(2,3,2,1) = 1;
- N = [15 5];
- m = [1 -1 -1];
- r = [1 1]; # reference station is station 1
- [U R Q X] = qncmmva(N,S,P,r,m);
- # FIXME: I replaced 14.3->14.323
- assert( U, [0.358 14.323 0; 0.156 0 4.707], 1e-3 );
- # FIXME: I replaced 14.3->14.323
- assert( X, [14.323 14.3230 0; .314 0 .314 ], 1e-3 );
- # FIXME: I replaced 14.3->14.323
- assert( Q, [.677 14.323 0; .293 0 4.707], 1e-3 );
- assert( R, [.047 1 15.0; .934 1 15.0], 1e-3 );
- assert( Q, R.*X, 1e-5 ); # Little's Law
-***** test
- C = 2; K = 3;
- S = [.01 .07 .10; ...
-      .05 .07 .10 ];
- P = zeros(C,K,C,K);
- P(1,1,1,2) = .7;
- P(1,1,1,3) = .2;
- P(1,1,2,1) = .1;
- P(2,1,2,2) = .3;
- P(2,1,2,3) = .5;
- P(2,1,1,1) = .2;
- P(1,2,1,1) = P(2,2,2,1) = 1;
- P(1,3,1,1) = P(2,3,2,1) = 1;
- N = [3 0];
- [U R Q X] = qncmmva(N, S, P);
- assert( R, [.015 .133 .163; .073 .133 .163], 1e-3 );
- assert( X, [12.609 8.826 2.522; 6.304 1.891 3.152], 1e-3 );
- assert( Q, [.185 1.175 .412; .462 .252 .515], 1e-3 );
- assert( U, [.126 .618 .252; .315 .132 .315], 1e-3 );
-***** test
- V = [ 1.00 0.45 0.50 0.00; ...
-       1.00 0.00 0.50 0.49 ];
- N = [3 2];
- S = [0.01 0.09 0.10 0.08; ...
-      0.05 0.09 0.10 0.08];
- [U R Q X] = qncmmva(N, S, V);
- assert( U, [ 0.1215 0.4921 0.6075 0.0000; ...
-              0.3433 0.0000 0.3433 0.2691 ], 1e-4 );
- assert( Q, [ 0.2131 0.7539 2.0328 0.0000; ...
-              0.5011 0.0000 1.1839 0.3149 ], 1e-4 );
- assert( R.*V, [0.0175 0.0620 0.1672 0.0000; ...
-                0.0729 0.0000 0.1724 0.0458 ], 1e-4 );
- assert( X, [12.1517 5.4682 6.0758 0.0000; ...
-              6.8669 0.0000 3.4334 3.3648 ], 1e-4 );
-***** test
- N = [1   1   1];
- S = [0.20000   0.02000   0.00000   0.00000;
-      0.20000   0.00000   0.02000   0.00000;
-      0.20000   0.00000   0.00000   0.02000];
- V = [1   1   0   0;
-      1   0   1   0;
-      1   0   0   1];
- [U R Q X] = qncmmva(N,S,V);
- assert( Q(1,2), Q(2,3), 1e-5);
- assert( Q(2,3), Q(3,4), 1e-5);
- assert( abs(max(Q(:,1)) - min(Q(:,1))) < 1e-5 );
-***** demo
- Ntot = 100; # total population size
- b = linspace(0.1,0.9,10); # fractions of class-1 requests
- S = [20 80 31 14 23 12; ...
-      90 30 33 20 14 7];
- V = ones(size(S));
- X1 = X1 = XX = zeros(size(b));
- R1 = R2 = RR = zeros(size(b));
- for i=1:length(b)
-   N = [fix(b(i)*Ntot) Ntot-fix(b(i)*Ntot)];
-   # printf("[%3d %3d]\n", N(1), N(2) );
-   [U R Q X] = qncmmva( N, S, V );
-   X1(i) = X(1,1) / V(1,1);
-   X2(i) = X(2,1) / V(2,1);
-   XX(i) = X1(i) + X2(i);
-   R1(i) = dot(R(1,:), V(1,:));
-   R2(i) = dot(R(2,:), V(2,:));
-   RR(i) = Ntot / XX(i);
- endfor
- subplot(2,1,1);
- plot(b, X1, ";Class 1;", "linewidth", 2, ...
-      b, X2, ";Class 2;", "linewidth", 2, ...
-      b, XX, ";System;", "linewidth", 2 );
- legend("location","south"); legend("boxoff");
- ylabel("Throughput");
- subplot(2,1,2);
- plot(b, R1, ";Class 1;", "linewidth", 2, ...
-      b, R2, ";Class 2;", "linewidth", 2, ...
-      b, RR, ";System;", "linewidth", 2 );
- legend("location","south"); legend("boxoff");
- xlabel("Population mix \\beta for Class 1");
- ylabel("Resp. Time");
-***** demo
- S = [1 0 .025; 0 15 .5];
- P = zeros(2,3,2,3);
- P(1,1,1,3) = P(1,3,1,1) = 1;
- P(2,2,2,3) = P(2,3,2,2) = 1;
- V = qncmvisits(P,[3 3]); # reference station is station 3
- N = [15 5];
- m = [-1 -1 1];
- [U R Q X] = qncmmva(N,S,V,m)
-***** demo
- C = 2; K = 3;
- S = [.01 .07 .10; ...
-      .05 .07 .10 ];
- P = zeros(C,K,C,K);
- P(1,1,1,2) = .7; P(1,1,1,3) = .2; P(1,1,2,1) = .1;
- P(2,1,2,2) = .3; P(2,1,2,3) = .5; P(2,1,1,1) = .2;
- P(1,2,1,1) = P(2,2,2,1) = 1;
- P(1,3,1,1) = P(2,3,2,1) = 1;
- N = [3 0];
- [U R Q X] = qncmmva(N, S, P)
-17 tests, 17 passed, 0 known failure, 0 skipped
-[inst/qsammm.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qsammm.m
-***** test
- [U R Q X] = qsammm( 73,[10,15,20,20,25] );
- assert( U, 0.81, 1e-2 );
- assert( Q, 6.5278, 1e-4 );
-1 test, 1 passed, 0 known failure, 0 skipped
-[inst/dtmcbd.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/dtmcbd.m
-***** test
- birth = [.5 .5 .3];
- death = [.6 .2 .3];
- fail("dtmcbd(birth,death)","must be");
-***** demo
- birth = [ .2 .3 .4 ];
- death = [ .1 .2 .3 ];
- P = dtmcbd( birth, death );
- disp(P)
-1 test, 1 passed, 0 known failure, 0 skipped
-[inst/qncsmvald.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsmvald.m
-***** test
- N = 0;
- S = [1 2; 3 4; 5 6; 7 8];
- V = [1 1 1 4];
- [U R Q X] = qncsmvald(N, S, V);
- assert( U, 0*V );
- assert( R, 0*V );
- assert( Q, 0*V );
- assert( X, 0*V );
-***** test
- # Exsample 3.42 p. 577 Jain
- V = [ 16 10 5 ];
- N = 20;
- S = [ 0.125 0.3 0.2 ];
- Sld = repmat( S', 1, N );
- Z = 4;
- [U1 R1 Q1 X1] = qncsmvald(N,Sld,V,Z);
- [U2 R2 Q2 X2] = qncsmva(N,S,V,ones(1,3),Z);
- assert( U1, U2, 1e-3 );
- assert( R1, R2, 1e-3 );
- assert( Q1, Q2, 1e-3 );
- assert( X1, X2, 1e-3 );
-***** test
- # Example 8.7 p. 349 Bolch et al.
- N = 3;
- Sld = 1 ./ [ 2 4 4; ...
-              1.667 1.667 1.667; ...
-              1.25 1.25 1.25; ...
-              1 2 3 ];
- V = [ 1 .5 .5 1 ];
- [U R Q X] = qncsmvald(N,Sld,V);
- assert( Q, [0.624 0.473 0.686 1.217], 1e-3 );
- assert( R, [0.512 0.776 1.127 1], 1e-3 );
- assert( all( U<=1) );
-***** test
- # Example 8.4 p. 333 Bolch et al.
- N = 3;
- S = [ .5 .6 .8 1 ];
- m = [2 1 1 N];
- Sld = zeros(3,N);
- Sld(1,:) = .5 ./ [1 2 2];
- Sld(2,:) = [.6 .6 .6];
- Sld(3,:) = [.8 .8 .8];
- Sld(4,:) = 1 ./ [1 2 3];
- V = [ 1 .5 .5 1 ];
- [U1 R1 Q1 X1] = qncsmvald(N,Sld,V);
- [U2 R2 Q2 X2] = qncsmva(N,S,V,m);
- ## Note that qncsmvald computes the utilization in a different
- ## way as qncsmva; in fact, qncsmva knows that service
- ## center i has m(i)>1 servers, but qncsmvald does not. Thus,
- ## utilizations for multiple-server nodes cannot be compared
- assert( U1([2,3]), U2([2,3]), 1e-3 );
- assert( R1, R2, 1e-3 );
- assert( Q1, Q2, 1e-3 );
- assert( X1, X2, 1e-3 );
-4 tests, 4 passed, 0 known failure, 0 skipped
-[inst/qnosbsb.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnosbsb.m
-***** test
- fail( "qnosbsb( 0.1, [] )", "vector" );
- fail( "qnosbsb( 0.1, [0 -1])", "nonnegative" );
- fail( "qnosbsb( 0, [1 2] )", "lambda" );
- fail( "qnosbsb( -1, [1 2])", "lambda" );
-***** test
- [Xl Xu Rl Ru] = qnosbsb(0.1,[1 2 3]);
- assert( Xl, 0 );
-2 tests, 2 passed, 0 known failure, 0 skipped
-[inst/qncsmva.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsmva.m
-***** test
- fail( "qncsmva()", "Invalid" );
- fail( "qncsmva( 10, [1 2], [1 2 3] )", "incompatible size" );
- fail( "qncsmva( 10, [-1 1], [1 1] )", "nonnegative" );
- fail( "qncsmva( 10.3, [-1 1], [1 1] )", "integer" );
- fail( "qncsmva( -0.3, [-1 1], [1 1] )", "nonnegative" );
-***** test
- qncsmva(1,1,1,1);
- qncsmva(1,1,1,-1);
- qncsmva(1,1,1,2);
- qncsmva(1,[1 1],[1 1],[-1 -1]);
- qncsmva(1,[1 1],[1 1],[1 1]);
- qncsmva(1,[1 1],[1 1],[2 2]);
-***** test
- N = 0;
- S = [1 2 3 4];
- V = [1 1 1 4];
- [U R Q X] = qncsmva(N, S, V);
- assert( U, 0*S );
- assert( R, 0*S );
- assert( Q, 0*S );
- assert( X, 0*S );
-***** test
- N = 0;
- S = [1 2 3 4];
- V = [1 1 1 4];
- m = [2 3 4 5];
- [U R Q X] = qncsmva(N, S, V, m);
- assert( U, 0*S );
- assert( R, 0*S );
- assert( Q, 0*S );
- assert( X, 0*S );
-***** test
- # Exsample 3.42 p. 577 Jain
- S = [ 0.125 0.3 0.2 ]';
- V = [ 16 10 5 ];
- N = 20;
- m = ones(1,3)';
- Z = 4;
- [U R Q X] = qncsmva(N,S,V,m,Z);
- assert( R, [ .373 4.854 .300 ], 1e-3 );
- assert( Q, [ 1.991 16.177 0.500 ], 1e-3 );
- assert( all( U>=0 ) );
- assert( all( U<=1 ) );
- assert( Q, R.*X, 1e-5 ); # Little's Law
-***** test
- # Exsample 3.42 p. 577 Jain
- S = [ 0.125 0.3 0.2 ];
- V = [ 16 10 5 ];
- N = 20;
- m = ones(1,3);
- Z = 4;
- [U R Q X] = qncsmva(N,S,V,m,Z);
- assert( R, [ .373 4.854 .300 ], 1e-3 );
- assert( Q, [ 1.991 16.177 0.500 ], 1e-3 );
- assert( all( U>=0 ) );
- assert( all( U<=1 ) );
- assert( Q, R.*X, 1e-5 ); # Little's Law
-***** test
- # Example 8.4 p. 333 Bolch et al.
- S = [ .5 .6 .8 1 ];
- N = 3;
- m = [2 1 1 -1];
- V = [ 1 .5 .5 1 ];
- [U R Q X] = qncsmva(N,S,V,m);
- assert( Q, [ 0.624 0.473 0.686 1.217 ], 1e-3 );
- assert( X, [ 1.218 0.609 0.609 1.218 ], 1e-3 );
- assert( all(U >= 0 ) );
- assert( all(U( m>0 ) <= 1 ) );
- assert( Q, R.*X, 1e-5 ); # Little's Law
-***** test
- # Example 8.3 p. 331 Bolch et al.
- # This is a single-class network, which however nothing else than
- # a special case of multiclass network
- S = [ 0.02 0.2 0.4 0.6 ];
- K = 6;
- V = [ 1 0.4 0.2 0.1 ];
- [U R Q X] = qncsmva(K, S, V);
- assert( U, [ 0.198 0.794 0.794 0.595 ], 1e-3 );
- assert( R, [ 0.025 0.570 1.140 1.244 ], 1e-3 );
- assert( Q, [ 0.244 2.261 2.261 1.234 ], 1e-3 );
- assert( X, [ 9.920 3.968 1.984 0.992 ], 1e-3 );
-***** test
- # Check bound analysis
- N = 10; # max population
- for n=1:N
-   S = [1 0.8 1.2 0.5];
-   V = [1 2 2 1];
-   [U R Q X] = qncsmva(n, S, V);
-   Xs = X(1)/V(1);
-   Rs = dot(R,V);
-   # Compare with balanced system bounds
-   [Xlbsb Xubsb Rlbsb Rubsb] = qncsbsb( n, S .* V );
-   assert( Xlbsb<=Xs );
-   assert( Xubsb>=Xs );
-   assert( Rlbsb<=Rs );
-   assert( Rubsb>=Rs );
-   # Compare with asymptotic bounds
-   [Xlab Xuab Rlab Ruab] = qncsaba( n, S .* V );
-   assert( Xlab<=Xs );
-   assert( Xuab>=Xs );
-   assert( Rlab<=Rs );
-   assert( Ruab>=Rs );
- endfor
-***** demo
- S = [ 0.125 0.3 0.2 ];
- V = [ 16 10 5 ];
- N = 20;
- m = ones(1,3);
- Z = 4;
- [U R Q X] = qncsmva(N,S,V,m,Z);
- X_s = X(1)/V(1); # System throughput
- R_s = dot(R,V); # System response time
- printf("\t    Util      Qlen     RespT      Tput\n");
- printf("\t--------  --------  --------  --------\n");
- for k=1:length(S)
-   printf("Dev%d\t%8.4f  %8.4f  %8.4f  %8.4f\n", k, U(k), Q(k), R(k), X(k) );
- endfor
- printf("\nSystem\t          %8.4f  %8.4f  %8.4f\n\n", N-X_s*Z, R_s, X_s );
-9 tests, 9 passed, 0 known failure, 0 skipped
-[inst/ctmcmtta.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/ctmcmtta.m
-***** test
- Q = [0 1 0; 1 0 1; 0 1 0 ]; Q -= diag( sum(Q,2) );
- fail( "ctmcmtta(Q,[1 0 0])", "no absorbing");
-***** test
- Q = [0 1 0; 1 0 1; 0 0 0; 0 0 0 ];
- fail( "ctmcmtta(Q,[1 0 0])", "square matrix");
-***** test
- Q = [0 1 0; 1 0 1; 0 0 0 ];
- fail( "ctmcmtta(Q,[1 0 0])", "infinitesimal");
-***** test
- Q = [ 0 0.1 0 0; ...
-       0.9 0 0.1 0; ...
-       0 0.9 0 0.1; ...
-       0 0 0 0 ];
- Q -= diag( sum(Q,2) );
- assert( ctmcmtta( Q,[0 0 0 1] ), 0 ); # state 4 is absorbing
-***** test
- Q = [-1 1; 0 0];
- assert( ctmcmtta( Q, [0 1] ), 0 ); # state 2 is absorbing
- assert( ctmcmtta( Q, [1 0] ), 1 ); # the result has been computed by hand
-***** demo
- mu = 0.01;
- death = [ 3 4 5 ] * mu;
- birth = 0*death;
- Q = ctmcbd(birth,death);
- t = ctmcmtta(Q,[0 0 0 1])
-***** demo
- N = 100;
- birth = death = ones(1,N-1); birth(1) = death(N-1) = 0;
- Q = diag(birth,1)+diag(death,-1); 
- Q -= diag(sum(Q,2));
- t = zeros(1,N/2);
- initial_state = 1:(N/2);
- for i=initial_state
-   p = zeros(1,N); p(i) = 1;
-   t(i) = ctmcmtta(Q,p);
- endfor
- plot(initial_state,t,"+");
- xlabel("Initial state");
- ylabel("MTTA");
-5 tests, 5 passed, 0 known failure, 0 skipped
-[inst/dtmcfpt.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/dtmcfpt.m
-***** test
- P = [1 1 1; 1 1 1];
- fail( "dtmcfpt(P)" );
-***** test
- P = dtmcbd([1 1 1], [0 0 0] );
- fail( "dtmcfpt(P)", "absorbing" );
-***** test
- P = [ 0.0 0.9 0.1; ...
-       0.1 0.0 0.9; ...
-       0.9 0.1 0.0 ];
- p = dtmc(P);
- M = dtmcfpt(P);
- assert( diag(M)', 1./p, 1e-8 );
-***** test
- P = [ 0 1 0 0 0; ...
-      .25 .0 .75 0 0; ...
-      0 .5 0 .5 0; ...
-      0 0 .75 0 .25; ...
-      0 0 0 1 0 ];
- M = dtmcfpt(P);
- assert( M, [16 1 2.6667 6.3333 21.3333; ...
-             15 4 1.6667 5.3333 20.3333; ...
-             18.6667 3.6667 2.6667 3.6667 18.6667; ...
-             20.3333 5.3333 1.6667 4 15; ...
-             21.3333 6.3333 2.6667 1 16 ], 1e-4 );
-***** test
- sz = 10;
- P = reshape( 1:sz^2, sz, sz );
- normP = repmat(sum(P,2),1,columns(P));
- P = P./normP;
- M = dtmcfpt(P);
- for i=1:rows(P)
-   for j=1:columns(P)
-     assert( M(i,j), 1 + dot(P(i,:), M(:,j)) - P(i,j)*M(j,j), 1e-8);
-   endfor
- endfor
-***** test
- P = zeros(9,9);
- P(1,[2 4]) = .5;
- P(2,[1 5 3]) = 1/3;
- P(3,[2 6]) = .5;
- P(4,[1 5 7]) = 1/3;
- P(5,[2 4 6 8]) = 1/4;
- P(6,[3 5 9]) = 1/3;
- P(7,[4 8]) = .5;
- P(8,[7 5 9]) = 1/3;
- P(9,[6 8]) = .5;
- M = dtmcfpt(P);
- assert( M(1:9 != 5,5)', [6 5 6 5 5 6 5 6], 100*eps );
-6 tests, 6 passed, 0 known failure, 0 skipped
 [inst/qsmm1.m]
 >>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qsmm1.m
 ***** test
@@ -4590,158 +4154,95 @@
  ylabel("p_k");
  title(sprintf("M/M/1 system, \\lambda = %g, \\mu = %g", lambda, mu));
 3 tests, 3 passed, 0 known failure, 0 skipped
-[inst/qnom.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnom.m
-***** test
- fail( "qnom([1 1], [.9; 1.0])", "exceeded at center 1");
- fail( "qnom([1 1], [0.9 .9; 0.9 1.0])", "exceeded at center 2");
- #qnom([1 1], [.9; 1.0],[],2); # should not fail, M/M/2-FCFS
- #qnom([1 1], [.9; 1.0],[],-1); # should not fail, -/G/1-PS
- fail( "qnom(1./[2 3], [1.9 1.9 0.9; 2.9 3.0 2.9])", "exceeded at center 2");
- #qnom(1./[2 3], [1 1.9 0.9; 0.3 3.0 1.5],[],[1 2 1]); # should not fail
-***** test
- V = [1 1; 1 1];
- S = [1 3; 2 4];
- lambda = [3/19 2/19];
- [U R Q X] = qnom(lambda, S, diag( lambda / sum(lambda) ) * V );
- assert( U(1,1), 3/19, 1e-6 );
- assert( U(2,1), 4/19, 1e-6 );
- assert( R(1,1), 19/12, 1e-6 );
- assert( R(1,2), 57/2, 1e-6 );
- assert( Q(1,1), .25, 1e-6 );
- assert( Q, R.*X, 1e-5 ); # Little's Law
-***** test
- # example p. 138 Zahorjan et al.
- V = [ 10 9; 5 4];
- S = [ 1/10 1/3; 2/5 1];
- lambda = [3/19 2/19];
- [U R Q X] = qnom(lambda, S, diag( lambda / sum(lambda) ) * V );
- assert( X(1,1), 1.58, 1e-2 );
- assert( U(1,1), .158, 1e-3 );
- assert( R(1,1), .158, 1e-3 ); # modified from the original example, as the reference above considers R as the residence time, not the response time
- assert( Q(1,1), .25, 1e-2 );
- assert( Q, R.*X, 1e-5 ); # Little's Law
+[inst/qsmmmk.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qsmmmk.m
 ***** test
- # example 7.7 p. 304 Bolch et al. Please note that the book uses the 
- # notation P(i,r,j,s) (i,j are service centers, r,s are job
- # classes) while the queueing package uses P(r,i,s,j)
- P = zeros(2,3,2,3);
- lambda = S = zeros(2,3);
- P(1,1,1,2) = 0.4;
- P(1,1,1,3) = 0.3;
- P(1,2,1,1) = 0.6;
- P(1,2,1,3) = 0.4;
- P(1,3,1,1) = 0.5;
- P(1,3,1,2) = 0.5;
- P(2,1,2,2) = 0.3;
- P(2,1,2,3) = 0.6;
- P(2,2,2,1) = 0.7;
- P(2,2,2,3) = 0.3;
- P(2,3,2,1) = 0.4;
- P(2,3,2,2) = 0.6;
- S(1,1) = 1/8;
- S(1,2) = 1/12;
- S(1,3) = 1/16;
- S(2,1) = 1/24;
- S(2,2) = 1/32;
- S(2,3) = 1/36;
- lambda(1,1) = lambda(2,1) = 1;
- V = qnomvisits(P,lambda);
- assert( V, [ 3.333 2.292 1.917; 10 8.049 8.415] ./ 2, 1e-3);
- [U R Q X] = qnom(sum(lambda,2), S, V);
- assert( sum(U,1), [0.833 0.442 0.354], 1e-3 );
- # Note: the value of K_22 (corresponding to Q(2,2)) reported in the book
- # is 0.5. However, hand computation using the exact same formulas
- # from the book produces a different value, 0.451
- assert( Q, [2.5 0.342 0.186; 2.5 0.451 0.362], 1e-3 );
+ lambda = mu = m = 1;
+ k = 10;
+ [U R Q X p0] = qsmmmk(lambda,mu,m,k);
+ assert( Q, k/2, 1e-7 );
+ assert( U, 1-p0, 1e-7 );
 ***** test
- P = zeros(2,2,2,2);
- P(1,1,1,2) = 0.8; P(1,2,1,1) = 1;
- P(2,1,2,2) = 0.9; P(2,2,2,1) = 1;
- S = zeros(2,2);
- S(1,1) = 1.5; S(1,2) = 1.2;
- S(2,1) = 0.8; S(2,2) = 2.5;
- lambda = zeros(2,2);
- lambda(1,1) = 1/20;
- lambda(2,1) = 1/30;
- [U1 R1 Q1 X1] = qnom(lambda, S, P); # qnom_cs
- [U2 R2 Q2 X2] = qnom(sum(lambda,2), S, qnomvisits(P,lambda)); # qnom_nocs
+ lambda = [1 0.8 2 9.2 0.01];
+ mu = lambda + 0.17;
+ k = 12;
+ [U1 R1 Q1 X1] = qsmm1k(lambda,mu,k);
+ [U2 R2 Q2 X2] = qsmmmk(lambda,mu,1,k);
  assert( U1, U2, 1e-5 );
  assert( R1, R2, 1e-5 );
  assert( Q1, Q2, 1e-5 );
  assert( X1, X2, 1e-5 );
-***** demo
- P = zeros(2,2,2,2);
- lambda = zeros(2,2);
- S = zeros(2,2);
- P(1,1,2,1) = P(1,2,2,1) = 0.2;
- P(1,1,2,2) = P(2,2,2,2) = 0.8;
- S(1,1) = S(1,2) = 0.1;
- S(2,1) = S(2,2) = 0.05;
- rr = 1:100;
- Xk = zeros(2,length(rr));
- for r=rr
-   lambda(1,1) = lambda(1,2) = 1/r;
-   [U R Q X] = qnom(lambda,S,P);
-   Xk(:,r) = sum(X,1)';
- endfor
- plot(rr,Xk(1,:),";Server 1;","linewidth",2, ...
-      rr,Xk(2,:),";Server 2;","linewidth",2);
- legend("boxoff");
- xlabel("Class 1 interarrival time");
- ylabel("Throughput");
-5 tests, 5 passed, 0 known failure, 0 skipped
-[inst/qncmpopmix.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncmpopmix.m
-***** demo
- N = [2 3];
- mix = qncmpopmix(3, N)
-***** test
- N = [-1 2 2];
- fail( "qncmpopmix(1, N)" );
-***** test
- N = [2 2 2];
- fail( "qncmpopmix(7, N)" );
-***** test
- N = [2 3 4];
- f = qncmpopmix(0, N );
- assert( f, [0 0 0] );
- f = qncmpopmix( 1, N );
- assert( f, [1 0 0; 0 1 0; 0 0 1] );
- f = qncmpopmix( 2, N );
- assert( f, [2 0 0; 1 1 0; 0 2 0; 1 0 1; 0 1 1; 0 0 2] );
- f = qncmpopmix( 3, N );
- assert( f, [2 1 0; 1 2 0; 0 3 0; 2 0 1; 1 1 1; 0 2 1; 1 0 2; 0 1 2; 0 0 3] );
+ #assert( [U1 R1 Q1 X1], [U2 R2 Q2 X2], 1e-5 );
 ***** test
- N = [2 1];
- f = qncmpopmix( 1, N );
- assert( f, [1 0; 0 1] );
- f = qncmpopmix( 2, N );
- assert( f,  [2 0; 1 1] );
-4 tests, 4 passed, 0 known failure, 0 skipped
-[inst/dtmcchkP.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/dtmcchkP.m
+ lambda = 0.9;
+ mu = 0.75;
+ k = 10;
+ [U1 R1 Q1 X1 p01] = qsmmmk(lambda,mu,1,k);
+ [U2 R2 Q2 X2 p02] = qsmm1k(lambda,mu,k);
+ assert( [U1 R1 Q1 X1 p01], [U2 R2 Q2 X2 p02], 1e-5 );
 ***** test
- [r err] = dtmcchkP( [1 1 1; 1 1 1] );
- assert( r, 0 );
- assert( index(err, "square") > 0 );
+ lambda = 0.8;
+ mu = 0.85;
+ m = 3;
+ k = 5;
+ [U1 R1 Q1 X1 p0] = qsmmmk( lambda, mu, m, k );
+ birth = lambda*ones(1,k);
+ death = [ mu*linspace(1,m,m) mu*m*ones(1,k-m) ];
+ q = ctmc(ctmcbd( birth, death ));
+ U2 = dot( q, min( 0:k, m )/m );
+ assert( U1, U2, 1e-4 );
+ Q2 = dot( [0:k], q );
+ assert( Q1, Q2, 1e-4 );
+ assert( p0, q(1), 1e-4 );
 ***** test
- [r err] = dtmcchkP( [1 0 0; 0 0.5 0; 0 0 0] );
- assert( r, 0 );
- assert( index(err, "stochastic") > 0 );
+ # This test comes from an example I found on the web 
+ lambda = 40;
+ mu = 30;
+ m = 3;
+ k = 7;
+ [U R Q X p0] = qsmmmk( lambda, mu, m, k );
+ assert( p0, 0.255037, 1e-6 );
+ assert( R, 0.036517, 1e-6 );
 ***** test
- P = [0 1; 1 0];
- assert( dtmcchkP(P), 2 );
+ # This test comes from an example I found on the web 
+ lambda = 50;
+ mu = 10;
+ m = 4;
+ k = 6;
+ [U R Q X p0 pk] = qsmmmk( lambda, mu, m, k );
+ assert( pk, 0.293543, 1e-6 );
 ***** test
- P = dtmcbd( linspace(0.1,0.4,10), linspace(0.4,0.1,10) );
- assert( dtmcchkP(P), rows(P) );
+ # This test comes from an example I found on the web 
+ lambda = 3;
+ mu = 2;
+ m = 2;
+ k = 5;
+ [U R Q X p0 pk] = qsmmmk( lambda, mu, m, k );
+ assert( p0, 0.179334, 1e-6 );
+ assert( pk, 0.085113, 1e-6 );
+ assert( Q, 2.00595, 1e-5 );
+ assert( R-1/mu, 0.230857, 1e-6 ); # waiting time in the queue
+***** demo
+ ## Given a M/M/m/K queue, compute the steady-state probability pn
+ ## of having n jobs in the systen.
+ lambda = 0.2;
+ mu = 0.25;
+ m = 5;
+ K = 20;
+ n = 0:10;
+ pn = qsmmmk(lambda, mu, m, K, n);
+ plot(n, pn, "-o", "linewidth", 2);
+ xlabel("N. of jobs (n)");
+ ylabel("P_n");
+ title(sprintf("M/M/%d/%d system, \\lambda = %g, \\mu = %g", m, K, lambda, mu));
+7 tests, 7 passed, 0 known failure, 0 skipped
+[inst/qsammm.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qsammm.m
 ***** test
- N = 1000;
- P = reshape( 1:N^2, N, N );
- P(1:N+1:end) = 0;
- P = P ./ repmat(sum(P,2),1,N);
- assert( dtmcchkP(P), N );
-5 tests, 5 passed, 0 known failure, 0 skipped
+ [U R Q X] = qsammm( 73,[10,15,20,20,25] );
+ assert( U, 0.81, 1e-2 );
+ assert( Q, 6.5278, 1e-4 );
+1 test, 1 passed, 0 known failure, 0 skipped
 [inst/qncmaba.m]
 >>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncmaba.m
 ***** test
@@ -4759,35 +4260,35 @@
     fail at line 130 column 7
     __test__ at line 3 column 2
     test at line 682 column 11
-    /tmp/tmp.cnygCL2JrV at line 166 column 31
+    /tmp/tmp.qG5gOSRbSf at line 118 column 31
 
 warning: the '.+' operator was deprecated in version 7 and will not be allowed in a future version of Octave; please use '+' instead; near line 124 of file '/build/reproducible-path/octave-queueing-1.2.7/debian/octave-queueing/usr/share/octave/packages/queueing-1.2.7/qncmaba.m'
 warning: called from
     fail at line 130 column 7
     __test__ at line 3 column 2
     test at line 682 column 11
-    /tmp/tmp.cnygCL2JrV at line 166 column 31
+    /tmp/tmp.qG5gOSRbSf at line 118 column 31
 
 warning: the '.+' operator was deprecated in version 7 and will not be allowed in a future version of Octave; please use '+' instead; near line 125 of file '/build/reproducible-path/octave-queueing-1.2.7/debian/octave-queueing/usr/share/octave/packages/queueing-1.2.7/qncmaba.m'
 warning: called from
     fail at line 130 column 7
     __test__ at line 3 column 2
     test at line 682 column 11
-    /tmp/tmp.cnygCL2JrV at line 166 column 31
+    /tmp/tmp.qG5gOSRbSf at line 118 column 31
 
 warning: the '.-' operator was deprecated in version 7 and will not be allowed in a future version of Octave; please use '-' instead; near line 126 of file '/build/reproducible-path/octave-queueing-1.2.7/debian/octave-queueing/usr/share/octave/packages/queueing-1.2.7/qncmaba.m'
 warning: called from
     fail at line 130 column 7
     __test__ at line 3 column 2
     test at line 682 column 11
-    /tmp/tmp.cnygCL2JrV at line 166 column 31
+    /tmp/tmp.qG5gOSRbSf at line 118 column 31
 
 warning: the '.-' operator was deprecated in version 7 and will not be allowed in a future version of Octave; please use '-' instead; near line 127 of file '/build/reproducible-path/octave-queueing-1.2.7/debian/octave-queueing/usr/share/octave/packages/queueing-1.2.7/qncmaba.m'
 warning: called from
     fail at line 130 column 7
     __test__ at line 3 column 2
     test at line 682 column 11
-    /tmp/tmp.cnygCL2JrV at line 166 column 31
+    /tmp/tmp.qG5gOSRbSf at line 118 column 31
 
 ***** test
  [Xl Xu Rl Ru] = qncmaba([0 0], [1 2 3; 1 2 3]);
@@ -4824,6 +4325,263 @@
  ylim([0, 700]);
  title("Class 2 response time"); legend("location", "northwest"); legend("boxoff");
 2 tests, 2 passed, 0 known failure, 0 skipped
+[inst/ctmcexps.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/ctmcexps.m
+***** test
+ Q = [-1 1; 1 -1];
+ L = ctmcexps(Q,10,[1 0]);
+ L = ctmcexps(Q,linspace(0,10,100),[1 0]);
+***** test
+ Q = ctmcbd( [1 2 3], [3 2 1] );
+ p0 = [1 0 0 0];
+ t = linspace(0,10,10);
+ L1 = L2 = zeros(length(t),4);
+ # compute L using the differential equation formulation
+ ff = @(x,t) (x(:)'*Q+p0);
+ fj = @(x,t) (Q);
+ L1 = lsode( {ff, fj}, zeros(size(p0)), t );
+ # compute L using ctmcexps (integral formulation)
+ for i=1:length(t)
+   L2(i,:) = ctmcexps(Q,t(i),p0);
+ endfor
+ assert( L1, L2, 1e-5);
+***** demo
+ lambda = 0.5;
+ N = 4;
+ b = lambda*[1:N-1];
+ d = zeros(size(b));
+ Q = ctmcbd(b,d);
+ t = linspace(0,10,100);
+ p0 = zeros(1,N); p0(1)=1;
+ L = zeros(length(t),N);
+ for i=1:length(t)
+   L(i,:) = ctmcexps(Q,t(i),p0);
+ endfor
+ plot( t, L(:,1), ";State 1;", "linewidth", 2, ...
+       t, L(:,2), ";State 2;", "linewidth", 2, ...
+       t, L(:,3), ";State 3;", "linewidth", 2, ...
+       t, L(:,4), ";State 4;", "linewidth", 2 );
+ legend("location","northwest"); legend("boxoff");
+ xlabel("Time");
+ ylabel("Expected sojourn time");
+***** demo
+ lambda = 0.5;
+ N = 4;
+ b = lambda*[1:N-1];
+ d = zeros(size(b));
+ Q = ctmcbd(b,d);
+ p0 = zeros(1,N); p0(1)=1;
+ L = ctmcexps(Q,p0);
+ disp(L);
+2 tests, 2 passed, 0 known failure, 0 skipped
+[inst/qncmmvabs.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncmmvabs.m
+***** test
+ S = [ 1 3 3; 2 4 3];
+ V = [ 1 1 3; 1 1 3];
+ N = [ 1 1 ];
+ m = [1 ; 1 ];
+ Z = [2 2 2];
+ fail( "qncmmvabs(N,S,V,m,Z)", "m must be" );
+ m = [1 ; 1 ; 1];
+ fail( "qncmmvabs(N,S,V,m,Z)", "Z must be" );
+warning: the '.+' operator was deprecated in version 7 and will not be allowed in a future version of Octave; please use '+' instead; near line 228 of file '/build/reproducible-path/octave-queueing-1.2.7/debian/octave-queueing/usr/share/octave/packages/queueing-1.2.7/qncmmvabs.m'
+warning: called from
+    fail at line 130 column 7
+    __test__ at line 8 column 2
+    test at line 682 column 11
+    /tmp/tmp.qG5gOSRbSf at line 134 column 31
+
+***** test
+ S = [ 1 3; 2 4];
+ V = [ 1 1; 1 1];
+ N = [ 1 1 ];
+ m = ones(1,2);
+ [U R Q X] = qncmmvabs(N,S,V,m);
+ assert( Q, [ .192 .808; .248 .752 ], 1e-3 );
+ Xc = ( X(:,1)./V(:,1) )';
+ assert( Xc, [ .154 .104 ], 1e-3 );
+ # Compute the (overall) class-c system response time
+ R_c = N ./ Xc;
+ assert( R_c, [ 6.508 9.614 ], 5e-3 );
+***** demo
+ S = [ 1, 1, 1, 1; 2, 1, 3, 1; 4, 2, 3, 3 ];
+ V = ones(3,4);
+ N = [10 5 1];
+ m = [1 0 1 1];
+ [U R Q X] = qncmmvabs(N,S,V,m);
+2 tests, 2 passed, 0 known failure, 0 skipped
+[inst/qnomvisits.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnomvisits.m
+***** test
+ fail( "qnomvisits( zeros(3,3,3), [1 1 1] )", "matrix");
+***** test
+ C = 2; K = 4;
+ P = zeros(C,K,C,K);
+ # class 1 routing
+ P(1,1,1,1) = .05;
+ P(1,1,1,2) = .45;
+ P(1,1,1,3) = .5;
+ P(1,2,1,1) = 0.1;
+ P(1,3,1,1) = 0.2;
+ # class 2 routing
+ P(2,1,2,1) = .01;
+ P(2,1,2,3) = .5;
+ P(2,1,2,4) = .49;
+ P(2,3,2,1) = 0.2;
+ P(2,4,2,1) = 0.16;
+ lambda = [0.1 0 0 0.1 ; 0 0 0.2 0.1];
+ lambda_sum = sum(lambda(:));
+ V = qnomvisits(P, lambda);
+ assert( all(V(:)>=0) );
+ for i=1:K
+   for c=1:C
+     assert(V(c,i), lambda(c,i) / lambda_sum + sum(sum(V .* P(:,:,c,i))), 1e-5);
+   endfor
+ endfor
+***** test
+ # example 7.7 p. 304 Bolch et al. 
+ # Note that the book uses a slightly different notation than
+ # what we use here. Specifically, the book defines the routing 
+ # probabilities as P(i,r,j,s) (i,j are service centers, r,s are job
+ # classes) while the queueing package uses P(r,i,s,j).
+ # A more serious problem arises in the definition of external arrivals.
+ # The computation of V(r,i) as given in the book (called e_ir
+ # in Eq 7.14) is performed in terms of P_{0, js}, defined as 
+ # "the probability in an open network that a job from outside the network
+ #  enters the jth node as a job of the sth class" (p. 267). This is
+ # compliant with eq. 7.12 where the external class r arrival rate at center
+ # i is computed as \lambda * P_{0,ir}. However, example 7.7 wrongly
+ # defines P_{0,11} = P_{0,12} = 1, instead of P_{0,11} = P_{0,12} = 0.5
+ # Therefore the resulting visit ratios they obtain must be divided by two.
+ P = zeros(2,3,2,3);
+ lambda = S = zeros(2,3);
+ P(1,1,1,2) = 0.4;
+ P(1,1,1,3) = 0.3;
+ P(1,2,1,1) = 0.6;
+ P(1,2,1,3) = 0.4;
+ P(1,3,1,1) = 0.5;
+ P(1,3,1,2) = 0.5;
+ P(2,1,2,2) = 0.3;
+ P(2,1,2,3) = 0.6;
+ P(2,2,2,1) = 0.7;
+ P(2,2,2,3) = 0.3;
+ P(2,3,2,1) = 0.4;
+ P(2,3,2,2) = 0.6;
+ lambda(1,1) = lambda(2,1) = 1;
+ V = qnomvisits(P,lambda);
+ assert( V, [ 3.333 2.292 1.917; 10 8.049 8.415] ./ 2, 1e-3);
+3 tests, 3 passed, 0 known failure, 0 skipped
+[inst/dtmcbd.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/dtmcbd.m
+***** test
+ birth = [.5 .5 .3];
+ death = [.6 .2 .3];
+ fail("dtmcbd(birth,death)","must be");
+***** demo
+ birth = [ .2 .3 .4 ];
+ death = [ .1 .2 .3 ];
+ P = dtmcbd( birth, death );
+ disp(P)
+1 test, 1 passed, 0 known failure, 0 skipped
+[inst/qncsgb.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsgb.m
+***** test
+ fail( "qncsgb( 1, [] )", "vector" );
+ fail( "qncsgb( 1, [0 -1])", "nonnegative" );
+ fail( "qncsgb( 0, [1 2] )", "> 0" );
+ fail( "qncsgb( -1, [1 2])", "nonnegative" );
+ fail( "qncsgb( 1, [1 2],1,[1 -1])", "single server" );
+warning: the '.-' operator was deprecated in version 7 and will not be allowed in a future version of Octave; please use '-' instead; near line 182 of file '/build/reproducible-path/octave-queueing-1.2.7/debian/octave-queueing/usr/share/octave/packages/queueing-1.2.7/qncsgb.m'
+warning: called from
+    fail at line 130 column 7
+    __test__ at line 3 column 2
+    test at line 682 column 11
+    /tmp/tmp.qG5gOSRbSf at line 158 column 31
+
+warning: the '.-' operator was deprecated in version 7 and will not be allowed in a future version of Octave; please use '-' instead; near line 184 of file '/build/reproducible-path/octave-queueing-1.2.7/debian/octave-queueing/usr/share/octave/packages/queueing-1.2.7/qncsgb.m'
+warning: called from
+    fail at line 130 column 7
+    __test__ at line 3 column 2
+    test at line 682 column 11
+    /tmp/tmp.qG5gOSRbSf at line 158 column 31
+
+***** # shared test function
+***** function test_gb( D, expected, Z=0 )
+ for i=1:rows(expected)
+   N = expected(i,1);
+   [X_lower X_upper Q_lower Q_upper] = qncsgb(N,D,1,1,Z);
+   X_exp_lower = expected(i,2);
+   X_exp_upper = expected(i,3);
+   assert( [N X_lower X_upper], [N X_exp_lower X_exp_upper], 1e-4 )
+ endfor
+***** xtest
+ # table IV
+ D = [ 0.1 0.1 0.09 0.08 ];
+ #            N  X_lower  X_upper
+ expected = [ 2  4.3040   4.3174; ...
+              5  6.6859   6.7524; ...
+              10 8.1521   8.2690; ...
+              20 9.0947   9.2431; ...
+              80 9.8233   9.8765 ];
+ test_gb(D, expected);
+!!!!! known failure
+ASSERT errors for:  assert ([N, X_lower, X_upper],[N, X_exp_lower, X_exp_upper],1e-4)
+
+  Location  |  Observed  |  Expected  |  Reason
+    (3)         4.3813       4.3174      Abs err 0.063891 exceeds tol 0.0001 by 0.06
+***** xtest
+ # table V
+ D = [ 0.1 0.1 0.09 0.08 ];
+ Z = 1;
+ #            N  X_lower  X_upper
+ expected = [ 2  1.4319   1.5195; ...
+              5  3.3432   3.5582; ...
+              10 5.7569   6.1410; ...
+              20 8.0856   8.6467; ...
+              80 9.7147   9.8594];
+ test_gb(D, expected, Z);
+!!!!! known failure
+ASSERT errors for:  assert ([N, X_lower, X_upper],[N, X_exp_lower, X_exp_upper],1e-4)
+
+  Location  |  Observed  |  Expected  |  Reason
+    (3)         1.5219       1.5195      Abs err 0.0023822 exceeds tol 0.0001 by 0.002
+***** test
+ P = [0 0.3 0.7; 1 0 0; 1 0 0];
+ S = [1 0.6 0.2];
+ m = ones(1,3);
+ V = qncsvisits(P);
+ Z = 2;
+ Nmax = 20;
+ tol = 1e-5; # compensate for numerical errors
+ ## Test case with Z>0
+ for n=1:Nmax
+   [X_gb_lower X_gb_upper NC NC Q_gb_lower Q_gb_upper] = qncsgb(n, S.*V, 1, 1, Z);
+   [U R Q X] = qnclosed( n, S, V, m, Z );
+   X_mva = X(1)/V(1);
+   assert( X_gb_lower <= X_mva+tol );
+   assert( X_gb_upper >= X_mva-tol );
+   assert( Q_gb_lower <= Q+tol ); # compensate for numerical errors
+   assert( Q_gb_upper >= Q-tol ); # compensate for numerical errors
+ endfor
+***** test
+ P = [0 0.3 0.7; 1 0 0; 1 0 0];
+ S = [1 0.6 0.2];
+ V = qncsvisits(P);
+ Nmax = 20;
+ tol = 1e-5; # compensate for numerical errors
+
+ ## Test case with Z=0
+ for n=1:Nmax
+   [X_gb_lower X_gb_upper NC NC Q_gb_lower Q_gb_upper] = qncsgb(n, S.*V);
+   [U R Q X] = qnclosed( n, S, V );
+   X_mva = X(1)/V(1);
+   assert( X_gb_lower <= X_mva+tol );
+   assert( X_gb_upper >= X_mva-tol );
+   assert( Q_gb_lower <= Q+tol );
+   assert( Q_gb_upper >= Q-tol );
+ endfor
+5 tests, 3 passed, 2 known failures, 0 skipped
 [inst/ctmc.m]
 >>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/ctmc.m
 ***** test
@@ -4973,40 +4731,190 @@
  MTBF = ctmcmtta(Q, p0) / hour;
  printf("System MTBF %.2f hours\n",MTBF);
 9 tests, 9 passed, 0 known failure, 0 skipped
-[inst/qnopen.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnopen.m
+[inst/qncsmvald.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsmvald.m
 ***** test
- # Example 34.1 p. 572
- lambda = 3;
- V = [16 7 8];
- S = [0.01 0.02 0.03];
- [U R Q X] = qnopen( lambda, S, V );
- assert( R, [0.0192 0.0345 0.107], 1e-2 );
- assert( U, [0.48 0.42 0.72], 1e-2 );
+ N = 0;
+ S = [1 2; 3 4; 5 6; 7 8];
+ V = [1 1 1 4];
+ [U R Q X] = qncsmvald(N, S, V);
+ assert( U, 0*V );
+ assert( R, 0*V );
+ assert( Q, 0*V );
+ assert( X, 0*V );
 ***** test
- V = [1 1; 1 1];
- S = [1 3; 2 4];
- lambda = [3/19 2/19];
- [U R Q] = qnopen(lambda, S, diag( lambda / sum(lambda) ) * V);
- assert( U(1,1), 3/19, 1e-6 );
- assert( U(2,1), 4/19, 1e-6 );
- assert( R(1,1), 19/12, 1e-6 );
- assert( R(1,2), 57/2, 1e-6 );
- assert( Q(1,1), .25, 1e-6 );
-2 tests, 2 passed, 0 known failure, 0 skipped
-[inst/ctmcbd.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/ctmcbd.m
+ # Exsample 3.42 p. 577 Jain
+ V = [ 16 10 5 ];
+ N = 20;
+ S = [ 0.125 0.3 0.2 ];
+ Sld = repmat( S', 1, N );
+ Z = 4;
+ [U1 R1 Q1 X1] = qncsmvald(N,Sld,V,Z);
+ [U2 R2 Q2 X2] = qncsmva(N,S,V,ones(1,3),Z);
+ assert( U1, U2, 1e-3 );
+ assert( R1, R2, 1e-3 );
+ assert( Q1, Q2, 1e-3 );
+ assert( X1, X2, 1e-3 );
 ***** test
- birth = [ 1 1 1 ];
- death = [ 2 2 2 ];
- Q = ctmcbd( birth, death );
- assert( ctmc(Q), [ 8/15 4/15 2/15 1/15 ], 1e-5 );
+ # Example 8.7 p. 349 Bolch et al.
+ N = 3;
+ Sld = 1 ./ [ 2 4 4; ...
+              1.667 1.667 1.667; ...
+              1.25 1.25 1.25; ...
+              1 2 3 ];
+ V = [ 1 .5 .5 1 ];
+ [U R Q X] = qncsmvald(N,Sld,V);
+ assert( Q, [0.624 0.473 0.686 1.217], 1e-3 );
+ assert( R, [0.512 0.776 1.127 1], 1e-3 );
+ assert( all( U<=1) );
+***** test
+ # Example 8.4 p. 333 Bolch et al.
+ N = 3;
+ S = [ .5 .6 .8 1 ];
+ m = [2 1 1 N];
+ Sld = zeros(3,N);
+ Sld(1,:) = .5 ./ [1 2 2];
+ Sld(2,:) = [.6 .6 .6];
+ Sld(3,:) = [.8 .8 .8];
+ Sld(4,:) = 1 ./ [1 2 3];
+ V = [ 1 .5 .5 1 ];
+ [U1 R1 Q1 X1] = qncsmvald(N,Sld,V);
+ [U2 R2 Q2 X2] = qncsmva(N,S,V,m);
+ ## Note that qncsmvald computes the utilization in a different
+ ## way as qncsmva; in fact, qncsmva knows that service
+ ## center i has m(i)>1 servers, but qncsmvald does not. Thus,
+ ## utilizations for multiple-server nodes cannot be compared
+ assert( U1([2,3]), U2([2,3]), 1e-3 );
+ assert( R1, R2, 1e-3 );
+ assert( Q1, Q2, 1e-3 );
+ assert( X1, X2, 1e-3 );
+4 tests, 4 passed, 0 known failure, 0 skipped
+[inst/erlangb.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/erlangb.m
+***** test
+ fail("erlangb(1, -1)", "positive");
+ fail("erlangb(-1, 1)", "positive");
+ fail("erlangb(1, 0)", "positive");
+ fail("erlangb(0, 1)", "positive");
+ fail("erlangb('foo',1)", "positive");
+ fail("erlangb(1,'bar')", "positive");
+ fail("erlangb([1 1],[1 1 1])","common size");
 1 test, 1 passed, 0 known failure, 0 skipped
-[inst/qnomaba.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnomaba.m
+[inst/qsmminf.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qsmminf.m
 ***** test
- fail( "qnomaba( [1 1], [1 1 1; 1 1 1; 1 1 1] )", "2 rows" );
+ fail( "qsmminf( [1 2], [1 2 3] )", "incompatible size");
+ fail( "qsmminf( [-1 -1], [1 1] )", ">0" );
+***** demo
+ ## Given a M/M/inf and M/M/m queue, compute the steady-state probability pk
+ ## of having k requests in the systen.
+ lambda = 5;
+ mu = 1.1;
+ m = 5;
+ k = 0:20;
+ pk_inf = qsmminf(lambda, mu, k);
+ pk_m = qsmmm(lambda, mu, 5, k);
+ plot(k, pk_inf, "-o;M/M/\\infty;", "linewidth", 2, ...
+      k, pk_m, "-x;M/M/5;", "linewidth", 2);
+ xlabel("N. of requests (k)");
+ ylabel("P_k");
+ title(sprintf("M/M/\\infty and M/M/%d systems, \\lambda = %g, \\mu = %g", m, lambda, mu));
 1 test, 1 passed, 0 known failure, 0 skipped
+[inst/qncspb.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncspb.m
+***** test
+ fail( "qncspb( 1, [] )", "vector" );
+ fail( "qncspb( 1, [0 -1])", "nonnegative" );
+ fail( "qncspb( 0, [1 2] )", "positive" );
+ fail( "qncspb( -1, [1 2])", "nonnegative" );
+ fail( "qncspb( 1, [1 2], [1,1], [2, 2])", "single server" );
+ fail( "qncspb( 1, [1 2], [1,1], [1, 1], -1)", "nonnegative" );
+***** # shared test function
+***** function test_pb( D, expected, Z=0 )
+ for i=1:rows(expected)
+   N = expected(i,1);
+   [X_lower X_upper] = qncspb(N,D,ones(size(D)),ones(size(D)),Z);
+   X_exp_lower = expected(i,2);
+   X_exp_upper = expected(i,3);
+   assert( [N X_lower X_upper], [N X_exp_lower X_exp_upper], 1e-4 )
+ endfor
+***** test
+ # table IV
+ D = [ 0.1 0.1 0.09 0.08 ];
+ #            N  X_lower  X_upper
+ expected = [ 2  4.3174   4.3174; ... 
+              5  6.6600   6.7297; ...
+              10 8.0219   8.2700; ...
+              20 8.8672   9.3387; ...
+              80 9.6736   10.000 ];
+ test_pb(D, expected);
+2 tests, 2 passed, 0 known failure, 0 skipped
+[inst/dtmcfpt.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/dtmcfpt.m
+***** test
+ P = [1 1 1; 1 1 1];
+ fail( "dtmcfpt(P)" );
+***** test
+ P = dtmcbd([1 1 1], [0 0 0] );
+ fail( "dtmcfpt(P)", "absorbing" );
+***** test
+ P = [ 0.0 0.9 0.1; ...
+       0.1 0.0 0.9; ...
+       0.9 0.1 0.0 ];
+ p = dtmc(P);
+ M = dtmcfpt(P);
+ assert( diag(M)', 1./p, 1e-8 );
+***** test
+ P = [ 0 1 0 0 0; ...
+      .25 .0 .75 0 0; ...
+      0 .5 0 .5 0; ...
+      0 0 .75 0 .25; ...
+      0 0 0 1 0 ];
+ M = dtmcfpt(P);
+ assert( M, [16 1 2.6667 6.3333 21.3333; ...
+             15 4 1.6667 5.3333 20.3333; ...
+             18.6667 3.6667 2.6667 3.6667 18.6667; ...
+             20.3333 5.3333 1.6667 4 15; ...
+             21.3333 6.3333 2.6667 1 16 ], 1e-4 );
+***** test
+ sz = 10;
+ P = reshape( 1:sz^2, sz, sz );
+ normP = repmat(sum(P,2),1,columns(P));
+ P = P./normP;
+ M = dtmcfpt(P);
+ for i=1:rows(P)
+   for j=1:columns(P)
+     assert( M(i,j), 1 + dot(P(i,:), M(:,j)) - P(i,j)*M(j,j), 1e-8);
+   endfor
+ endfor
+***** test
+ P = zeros(9,9);
+ P(1,[2 4]) = .5;
+ P(2,[1 5 3]) = 1/3;
+ P(3,[2 6]) = .5;
+ P(4,[1 5 7]) = 1/3;
+ P(5,[2 4 6 8]) = 1/4;
+ P(6,[3 5 9]) = 1/3;
+ P(7,[4 8]) = .5;
+ P(8,[7 5 9]) = 1/3;
+ P(9,[6 8]) = .5;
+ M = dtmcfpt(P);
+ assert( M(1:9 != 5,5)', [6 5 6 5 5 6 5 6], 100*eps );
+6 tests, 6 passed, 0 known failure, 0 skipped
+[inst/qnosaba.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnosaba.m
+***** test
+ fail( "qnosaba( 0.1, [] )", "vector" );
+ fail( "qnosaba( 0.1, [0 -1])", "nonnegative" );
+ fail( "qnosaba( 0, [1 2] )", "lambda" );
+ fail( "qnosaba( -1, [1 2])", "lambda" );
+ fail( "qnosaba( 1, [1 2 3], [1 2] )", "incompatible size");
+ fail( "qnosaba( 1, [1 2 3], [-1 2 3] )", "nonnegative");
+***** test
+ [Xl Xu Rl Ru] = qnosaba( 1, [1 1] );
+ assert( Xl, 0 );
+ assert( Ru, +inf );
+2 tests, 2 passed, 0 known failure, 0 skipped
 [inst/qsmm1k.m]
 >>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qsmm1k.m
 ***** test
@@ -5049,42 +4957,35 @@
  ylabel("p_n");
  title(sprintf("M/M/1/%d system, \\lambda = %g, \\mu = %g", K, lambda, mu));
 3 tests, 3 passed, 0 known failure, 0 skipped
-[inst/qncmmvabs.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncmmvabs.m
+[inst/qncsaba.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsaba.m
 ***** test
- S = [ 1 3 3; 2 4 3];
- V = [ 1 1 3; 1 1 3];
- N = [ 1 1 ];
- m = [1 ; 1 ];
- Z = [2 2 2];
- fail( "qncmmvabs(N,S,V,m,Z)", "m must be" );
- m = [1 ; 1 ; 1];
- fail( "qncmmvabs(N,S,V,m,Z)", "Z must be" );
-warning: the '.+' operator was deprecated in version 7 and will not be allowed in a future version of Octave; please use '+' instead; near line 228 of file '/build/reproducible-path/octave-queueing-1.2.7/debian/octave-queueing/usr/share/octave/packages/queueing-1.2.7/qncmmvabs.m'
-warning: called from
-    fail at line 130 column 7
-    __test__ at line 8 column 2
-    test at line 682 column 11
-    /tmp/tmp.cnygCL2JrV at line 214 column 31
-
+ fail("qncsaba(-1,0)", "N must be");
+ fail("qncsaba(1,[])", "nonempty");
+ fail("qncsaba(1,[-1 2])", "nonnegative");
+ fail("qncsaba(1,[1 2],[1 2 3])", "incompatible size");
+ fail("qncsaba(1,[1 2 3],[1 2 -1])", "nonnegative");
+ fail("qncsaba(1,[1 2 3],[1 2 3],[1 2])", "incompatible size");
+ fail("qncsaba(1,[1 2 3],[1 2 3],[1 2 1])", "not supported");
+ fail("qncsaba(1,[1 2 3],[1 2 3],[1 1 1],-1)", "nonnegative");
+ fail("qncsaba(1,[1 2 3],[1 2 3],1,[0 0])", "scalar");
 ***** test
- S = [ 1 3; 2 4];
- V = [ 1 1; 1 1];
- N = [ 1 1 ];
- m = ones(1,2);
- [U R Q X] = qncmmvabs(N,S,V,m);
- assert( Q, [ .192 .808; .248 .752 ], 1e-3 );
- Xc = ( X(:,1)./V(:,1) )';
- assert( Xc, [ .154 .104 ], 1e-3 );
- # Compute the (overall) class-c system response time
- R_c = N ./ Xc;
- assert( R_c, [ 6.508 9.614 ], 5e-3 );
-***** demo
- S = [ 1, 1, 1, 1; 2, 1, 3, 1; 4, 2, 3, 3 ];
- V = ones(3,4);
- N = [10 5 1];
- m = [1 0 1 1];
- [U R Q X] = qncmmvabs(N,S,V,m);
+ N = 20;
+ S = [ 4.6*2 8 ];
+ Z = 120;
+ [X_l X_u R_l R_u] = qncsaba(N, S, ones(size(S)), ones(size(S)), Z);
+ assert( [X_u R_l], [0.109 64], 1e-3 );
+2 tests, 2 passed, 0 known failure, 0 skipped
+[inst/qnosbsb.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnosbsb.m
+***** test
+ fail( "qnosbsb( 0.1, [] )", "vector" );
+ fail( "qnosbsb( 0.1, [0 -1])", "nonnegative" );
+ fail( "qnosbsb( 0, [1 2] )", "lambda" );
+ fail( "qnosbsb( -1, [1 2])", "lambda" );
+***** test
+ [Xl Xu Rl Ru] = qnosbsb(0.1,[1 2 3]);
+ assert( Xl, 0 );
 2 tests, 2 passed, 0 known failure, 0 skipped
 [inst/ctmctaexps.m]
 >>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/ctmctaexps.m
@@ -5155,6 +5056,489 @@
  axis(ax);
  legend("boxoff");
 1 test, 1 passed, 0 known failure, 0 skipped
+[inst/dtmcexps.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/dtmcexps.m
+***** test
+ P = dtmcbd([1 1 1 1], [0 0 0 0]);
+ L = dtmcexps(P,[1 0 0 0 0]);
+ t = dtmcmtta(P,[1 0 0 0 0]);
+ assert( L, [1 1 1 1 0] );
+ assert( sum(L), t );
+***** test
+ P = dtmcbd(linspace(0.1,0.4,5),linspace(0.4,0.1,5));
+ p0 = [1 0 0 0 0 0];
+ L = dtmcexps(P,0,p0);
+ assert( L, p0 );
+2 tests, 2 passed, 0 known failure, 0 skipped
+[inst/qnmix.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnmix.m
+***** test
+ lambda = [1 0 0];
+ N = [1 1 1];
+ S = V = [1 1 1; 1 1 1; 1 1 1];
+ fail( "qnmix( lambda, N, S, V)", "same time");
+ N = [0 0 1];
+ fail( "qnmix( lambda, N, S, V)", "open nor closed" );
+ N = [0 1 2];
+ m = [ 1 1 2 ];
+ fail( "qnmix( lambda, N, S, V, m)", "single-server and delay" );
+ S = V = [1 1 1; 1 1 1];
+ fail( "qnmix( lambda, N, S, V)", "rows" );
+***** test
+ # Example p. 148 Zahorjan et al.
+ lambda = [1 1/2 0 0];
+ N = [0 0 1 1];
+ V = [1 1; 1 1; 1 1; 1 1];
+ S = [1/4 1/6; 1/2 1; 1/2 1; 1 4/3];
+ [U R Q X] = qnmix(lambda, N, S, V );
+ assert( Q(3,1), 4/19, 1e-4 );
+ assert( Q(3,2), 15/19, 1e-4 );
+ assert( Q(4,1), 5/19, 1e-4 );
+ assert( Q(4,2), 14/19, 1e-4 );
+ assert( Q, R.*X, 1e-5 ); # Little's Law
+***** test
+ # Example 8.6 p. 345 Bolch et al.
+ lambda = [0.5 0.25 0 0];
+ N = [0 0 1 1];
+ V = [2 1; 2.5 1.5; 1 0.5; 1 0.4];
+ S = [0.4 0.6; 0.8 1.6; 0.3 0.5; 0.5 0.8];
+ [U R Q X] = qnmix( lambda, N, S, V );
+ assert( U([1 2],:), [0.4 0.3; 0.5 0.6], 1e-3 );
+ assert( R([3 4],:), [4.829 6.951; 7.727 11.636], 1e-3 );
+ assert( Q([3 4],:), [0.582 0.418; 0.624 0.376], 1e-3 );
+ assert( Q([1 2],:), [8.822 5.383; 11.028 10.766], 1e-3 ); #FIXME
+ assert( R([1 2],:), [8.822 10.766; 17.645 28.710], 1e-3 );
+ assert( X(3,1)/V(3,1), 0.120, 1e-3 );
+ assert( X(4,1)/V(4,1), 0.081, 1e-3 );
+ assert( Q, R.*X, 1e-5 ); # Little's Law
+***** test
+ ## example figure 10 p. 26 Schwetman, "Implementing the Mean Value
+ ## Analysis for the Solution of Queueing Network Models", Technical
+ ## Report CSD-TR-355, feb 15, 1982, Purdue University. 
+ S = [.25 0; .25 .10];
+ V = [1 0; 1 1];
+ lambda = [1 0];
+ N = [0 3];
+ [U R Q X] = qnmix( lambda, N, S, V );
+ assert( U(1,1), .25, 1e-3 );
+ assert( X(1,1), 1.0, 1e-3 );
+ assert( [R(1,1) R(2,1) R(2,2)], [1.201 0.885 0.135], 1e-3 );
+ assert( Q, R.*X, 1e-5 ); # Little's Law
+4 tests, 4 passed, 0 known failure, 0 skipped
+[inst/dtmc.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/dtmc.m
+***** test
+ P = [0.75 0.25; 0.5 0.5];
+ p = dtmc(P);
+ assert( p*P, p, 1e-5 );
+ assert( p, [0.6666 0.3333], 1e-4 );
+***** test
+ #Example 2.11 p. 44 Bolch et al.
+ P = [0.5 0.5; 0.5 0.5];
+ p = dtmc(P);
+ assert( p, [0.5 0.5], 1e-3 );
+***** test
+ fail("dtmc( [1 1 1; 1 1 1] )", "square");
+***** test
+ a = 0.2;
+ b = 0.8;
+ P = [1-a a; b 1-b];
+ plim = dtmc(P);
+ p = dtmc(P, 100, [1 0]);
+ assert( plim, p, 1e-5 );
+***** test
+ P = [0 1 0 0 0; ...
+      .25 0 .75 0 0; ...
+      0 .5 0 .5 0; ...
+      0 0 .75 0 .25; ...
+      0 0 0 1 0 ];
+ p = dtmc(P);
+ assert( p, [.0625 .25 .375 .25 .0625], 10*eps );
+***** test
+ P = zeros(9,9);
+ P(1,[2 4]) = 1/2;
+ P(2,[1 5 3]) = 1/3;
+ P(3,[2 6]) = 1/2;
+ P(4,[1 5 7]) = 1/3;
+ P(5,[2 4 6 8]) = 1/4;
+ P(6,[3 5 9]) = 1/3;
+ P(7,[4 8]) = 1/2;
+ P(8,[7 5 9]) = 1/3;
+ P(9,[6 8]) = 1/2;
+ p = dtmc(P);
+ assert( p, [1/12 1/8 1/12 1/8 1/6 1/8 1/12 1/8 1/12], 10*eps );
+***** demo
+ P = zeros(9,9);
+ P(1,[2 4]    ) = 1/2;
+ P(2,[1 5 3]  ) = 1/3;
+ P(3,[2 6]    ) = 1/2;
+ P(4,[1 5 7]  ) = 1/3;
+ P(5,[2 4 6 8]) = 1/4;
+ P(6,[3 5 9]  ) = 1/3;
+ P(7,[4 8]    ) = 1/2;
+ P(8,[7 5 9]  ) = 1/3;
+ P(9,[6 8]    ) = 1/2;
+ p = dtmc(P);
+ disp(p)
+***** demo
+ a = 0.2;
+ b = 0.15;
+ P = [ 1-a a; b 1-b];
+ T = 0:14;
+ pp = zeros(2,length(T));
+ for i=1:length(T)
+   pp(:,i) = dtmc(P,T(i),[1 0]);
+ endfor
+ ss = dtmc(P); # compute steady state probabilities
+ plot( T, pp(1,:), "b+;p_0(t);", "linewidth", 2, ...
+       T, ss(1)*ones(size(T)), "b;Steady State;", ...
+       T, pp(2,:), "r+;p_1(t);", "linewidth", 2, ...
+       T, ss(2)*ones(size(T)), "r;Steady State;" );
+ xlabel("Time Step");
+ legend("boxoff");
+6 tests, 6 passed, 0 known failure, 0 skipped
+[inst/qnomaba.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnomaba.m
+***** test
+ fail( "qnomaba( [1 1], [1 1 1; 1 1 1; 1 1 1] )", "2 rows" );
+1 test, 1 passed, 0 known failure, 0 skipped
+[inst/qncmpopmix.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncmpopmix.m
+***** demo
+ N = [2 3];
+ mix = qncmpopmix(3, N)
+***** test
+ N = [-1 2 2];
+ fail( "qncmpopmix(1, N)" );
+***** test
+ N = [2 2 2];
+ fail( "qncmpopmix(7, N)" );
+***** test
+ N = [2 3 4];
+ f = qncmpopmix(0, N );
+ assert( f, [0 0 0] );
+ f = qncmpopmix( 1, N );
+ assert( f, [1 0 0; 0 1 0; 0 0 1] );
+ f = qncmpopmix( 2, N );
+ assert( f, [2 0 0; 1 1 0; 0 2 0; 1 0 1; 0 1 1; 0 0 2] );
+ f = qncmpopmix( 3, N );
+ assert( f, [2 1 0; 1 2 0; 0 3 0; 2 0 1; 1 1 1; 0 2 1; 1 0 2; 0 1 2; 0 0 3] );
+***** test
+ N = [2 1];
+ f = qncmpopmix( 1, N );
+ assert( f, [1 0; 0 1] );
+ f = qncmpopmix( 2, N );
+ assert( f,  [2 0; 1 1] );
+4 tests, 4 passed, 0 known failure, 0 skipped
+[inst/ctmcmtta.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/ctmcmtta.m
+***** test
+ Q = [0 1 0; 1 0 1; 0 1 0 ]; Q -= diag( sum(Q,2) );
+ fail( "ctmcmtta(Q,[1 0 0])", "no absorbing");
+***** test
+ Q = [0 1 0; 1 0 1; 0 0 0; 0 0 0 ];
+ fail( "ctmcmtta(Q,[1 0 0])", "square matrix");
+***** test
+ Q = [0 1 0; 1 0 1; 0 0 0 ];
+ fail( "ctmcmtta(Q,[1 0 0])", "infinitesimal");
+***** test
+ Q = [ 0 0.1 0 0; ...
+       0.9 0 0.1 0; ...
+       0 0.9 0 0.1; ...
+       0 0 0 0 ];
+ Q -= diag( sum(Q,2) );
+ assert( ctmcmtta( Q,[0 0 0 1] ), 0 ); # state 4 is absorbing
+***** test
+ Q = [-1 1; 0 0];
+ assert( ctmcmtta( Q, [0 1] ), 0 ); # state 2 is absorbing
+ assert( ctmcmtta( Q, [1 0] ), 1 ); # the result has been computed by hand
+***** demo
+ mu = 0.01;
+ death = [ 3 4 5 ] * mu;
+ birth = 0*death;
+ Q = ctmcbd(birth,death);
+ t = ctmcmtta(Q,[0 0 0 1])
+***** demo
+ N = 100;
+ birth = death = ones(1,N-1); birth(1) = death(N-1) = 0;
+ Q = diag(birth,1)+diag(death,-1); 
+ Q -= diag(sum(Q,2));
+ t = zeros(1,N/2);
+ initial_state = 1:(N/2);
+ for i=initial_state
+   p = zeros(1,N); p(i) = 1;
+   t(i) = ctmcmtta(Q,p);
+ endfor
+ plot(initial_state,t,"+");
+ xlabel("Initial state");
+ ylabel("MTTA");
+5 tests, 5 passed, 0 known failure, 0 skipped
+[inst/qncsmva.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsmva.m
+***** test
+ fail( "qncsmva()", "Invalid" );
+ fail( "qncsmva( 10, [1 2], [1 2 3] )", "incompatible size" );
+ fail( "qncsmva( 10, [-1 1], [1 1] )", "nonnegative" );
+ fail( "qncsmva( 10.3, [-1 1], [1 1] )", "integer" );
+ fail( "qncsmva( -0.3, [-1 1], [1 1] )", "nonnegative" );
+***** test
+ qncsmva(1,1,1,1);
+ qncsmva(1,1,1,-1);
+ qncsmva(1,1,1,2);
+ qncsmva(1,[1 1],[1 1],[-1 -1]);
+ qncsmva(1,[1 1],[1 1],[1 1]);
+ qncsmva(1,[1 1],[1 1],[2 2]);
+***** test
+ N = 0;
+ S = [1 2 3 4];
+ V = [1 1 1 4];
+ [U R Q X] = qncsmva(N, S, V);
+ assert( U, 0*S );
+ assert( R, 0*S );
+ assert( Q, 0*S );
+ assert( X, 0*S );
+***** test
+ N = 0;
+ S = [1 2 3 4];
+ V = [1 1 1 4];
+ m = [2 3 4 5];
+ [U R Q X] = qncsmva(N, S, V, m);
+ assert( U, 0*S );
+ assert( R, 0*S );
+ assert( Q, 0*S );
+ assert( X, 0*S );
+***** test
+ # Exsample 3.42 p. 577 Jain
+ S = [ 0.125 0.3 0.2 ]';
+ V = [ 16 10 5 ];
+ N = 20;
+ m = ones(1,3)';
+ Z = 4;
+ [U R Q X] = qncsmva(N,S,V,m,Z);
+ assert( R, [ .373 4.854 .300 ], 1e-3 );
+ assert( Q, [ 1.991 16.177 0.500 ], 1e-3 );
+ assert( all( U>=0 ) );
+ assert( all( U<=1 ) );
+ assert( Q, R.*X, 1e-5 ); # Little's Law
+***** test
+ # Exsample 3.42 p. 577 Jain
+ S = [ 0.125 0.3 0.2 ];
+ V = [ 16 10 5 ];
+ N = 20;
+ m = ones(1,3);
+ Z = 4;
+ [U R Q X] = qncsmva(N,S,V,m,Z);
+ assert( R, [ .373 4.854 .300 ], 1e-3 );
+ assert( Q, [ 1.991 16.177 0.500 ], 1e-3 );
+ assert( all( U>=0 ) );
+ assert( all( U<=1 ) );
+ assert( Q, R.*X, 1e-5 ); # Little's Law
+***** test
+ # Example 8.4 p. 333 Bolch et al.
+ S = [ .5 .6 .8 1 ];
+ N = 3;
+ m = [2 1 1 -1];
+ V = [ 1 .5 .5 1 ];
+ [U R Q X] = qncsmva(N,S,V,m);
+ assert( Q, [ 0.624 0.473 0.686 1.217 ], 1e-3 );
+ assert( X, [ 1.218 0.609 0.609 1.218 ], 1e-3 );
+ assert( all(U >= 0 ) );
+ assert( all(U( m>0 ) <= 1 ) );
+ assert( Q, R.*X, 1e-5 ); # Little's Law
+***** test
+ # Example 8.3 p. 331 Bolch et al.
+ # This is a single-class network, which however nothing else than
+ # a special case of multiclass network
+ S = [ 0.02 0.2 0.4 0.6 ];
+ K = 6;
+ V = [ 1 0.4 0.2 0.1 ];
+ [U R Q X] = qncsmva(K, S, V);
+ assert( U, [ 0.198 0.794 0.794 0.595 ], 1e-3 );
+ assert( R, [ 0.025 0.570 1.140 1.244 ], 1e-3 );
+ assert( Q, [ 0.244 2.261 2.261 1.234 ], 1e-3 );
+ assert( X, [ 9.920 3.968 1.984 0.992 ], 1e-3 );
+***** test
+ # Check bound analysis
+ N = 10; # max population
+ for n=1:N
+   S = [1 0.8 1.2 0.5];
+   V = [1 2 2 1];
+   [U R Q X] = qncsmva(n, S, V);
+   Xs = X(1)/V(1);
+   Rs = dot(R,V);
+   # Compare with balanced system bounds
+   [Xlbsb Xubsb Rlbsb Rubsb] = qncsbsb( n, S .* V );
+   assert( Xlbsb<=Xs );
+   assert( Xubsb>=Xs );
+   assert( Rlbsb<=Rs );
+   assert( Rubsb>=Rs );
+   # Compare with asymptotic bounds
+   [Xlab Xuab Rlab Ruab] = qncsaba( n, S .* V );
+   assert( Xlab<=Xs );
+   assert( Xuab>=Xs );
+   assert( Rlab<=Rs );
+   assert( Ruab>=Rs );
+ endfor
+***** demo
+ S = [ 0.125 0.3 0.2 ];
+ V = [ 16 10 5 ];
+ N = 20;
+ m = ones(1,3);
+ Z = 4;
+ [U R Q X] = qncsmva(N,S,V,m,Z);
+ X_s = X(1)/V(1); # System throughput
+ R_s = dot(R,V); # System response time
+ printf("\t    Util      Qlen     RespT      Tput\n");
+ printf("\t--------  --------  --------  --------\n");
+ for k=1:length(S)
+   printf("Dev%d\t%8.4f  %8.4f  %8.4f  %8.4f\n", k, U(k), Q(k), R(k), X(k) );
+ endfor
+ printf("\nSystem\t          %8.4f  %8.4f  %8.4f\n\n", N-X_s*Z, R_s, X_s );
+9 tests, 9 passed, 0 known failure, 0 skipped
+[inst/qncmvisits.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncmvisits.m
+***** test
+
+ ## Closed, multiclass network
+
+ C = 2; K = 3; 
+ P = zeros(C,K,C,K);
+ P(1,1,1,2) = 1;
+ P(1,2,1,1) = 1;
+ P(2,1,2,3) = 1;
+ P(2,3,2,1) = 1;
+ V = qncmvisits(P);
+ for c=1:C
+   for k=1:K
+     assert(V(c,k), sum(sum(V .* P(:,:,c,k))), 1e-5);
+   endfor
+ endfor
+***** test
+
+ ## Test multiclass network. Example from Schwetman (figure 7, page 9 of
+ ## http://docs.lib.purdue.edu/cstech/259/
+ ## "Testing network-of-queues software, technical report CSD-TR 330,
+ ## Purdue University).
+
+ C = 2; K = 4;
+ P = zeros(C,K,C,K);
+ # class 1 routing
+ P(1,1,1,1) = .05;
+ P(1,1,1,2) = .45;
+ P(1,1,1,3) = .5;
+ P(1,2,1,1) = 1;
+ P(1,3,1,1) = 1;
+ # class 2 routing
+ P(2,1,2,1) = .01;
+ P(2,1,2,3) = .5;
+ P(2,1,2,4) = .49;
+ P(2,3,2,1) = 1;
+ P(2,4,2,1) = 1;
+ V = qncmvisits(P);
+ for c=1:C
+   for i=1:K
+     assert(V(c,i), sum(sum(V .* P(:,:,c,i))), 1e-5);
+   endfor
+ endfor
+***** test
+
+ ## Network with class switching.
+ ## This is the example in figure 9 of
+ ## Schwetman, "Implementing the Mean Value Analysis
+ ## Algorithm fort the solution of Queueing Network Models", Technical
+ ## Report CSD-TR-355, Department of Computer Science, Purdue Univrsity,
+ ## Feb 15, 1982, http://docs.lib.purdue.edu/cstech/286/
+
+ C = 2; K = 3;
+ S = [.01 .07 .10; ...
+      .05 0.7 .10 ];
+ P = zeros(C,K,C,K);
+ P(1,1,1,2) = .7;
+ P(1,1,1,3) = .2;
+ P(1,1,2,1) = .1;
+ P(2,1,2,2) = .3;
+ P(2,1,2,3) = .5;
+ P(2,1,1,1) = .2;
+ P(1,2,1,1) = P(1,3,1,1) = 1;
+ P(2,2,2,1) = P(2,3,2,1) = 1;
+ N = [3 0];
+ V = qncmvisits(P);
+ VV = [10 7 2; 5 1.5 2.5]; # result given in Schwetman; our function computes something different, but that's ok since visit counts are actually ratios
+ assert( V ./ repmat(V(:,1),1,K), VV ./ repmat(VV(:,1),1,K), 1e-5 );
+***** test
+
+ ## two disjoint classes: must produce two disjoing chains
+
+ C = 2; K = 3;
+ P = zeros(C,K,C,K);
+ P(1,1,1,2) = 1;
+ P(1,2,1,1) = 1;
+ P(2,1,2,3) = 1;
+ P(2,3,2,1) = 1;
+ [nc r] = qncmvisits(P);
+ assert( r(1) != r(2) );
+***** test
+
+ ## two classes, one chain
+
+ C = 2; K = 3;
+ P = zeros(C,K,C,K);
+ P(1,1,1,2) = .5;
+ P(1,2,2,1) = 1;
+ P(2,1,2,3) = .5;
+ P(2,3,1,1) = 1;
+ [nc r] = qncmvisits(P);
+ assert( r(1) == r(2) );
+***** test
+ 
+ ## a "Moebius strip". Note that this configuration is invalid, and
+ ## therefore our algorithm must raise an error. This is because this
+ ## network has two chains, but both chains contain both classes
+
+ C = 2; K = 2;
+ P = zeros(C,K,C,K);
+ P(1,1,2,2) = 1;
+ P(2,2,1,1) = 1;
+ P(2,1,1,2) = 1;
+ P(1,2,2,1) = 1;
+ fail( "qncmvisits(P)", "different");
+***** test
+
+ ## Network with two classes representing independent chains.
+ ## This is example in figure 8 of
+ ## Schwetman, "Implementing the Mean Value Analysis
+ ## Algorithm fort the solution of Queueing Network Models", Technical
+ ## Report CSD-TR-355, Department of Computer Science, Purdue Univrsity,
+ ## Feb 15, 1982, http://docs.lib.purdue.edu/cstech/286/
+ 
+ C = 2; K = 2;
+ P = zeros(C,K,C,K);
+ P(1,1,1,3) = P(1,3,1,1) = 1;
+ P(2,2,2,3) = P(2,3,2,2) = 1;
+ V = qncmvisits(P,[1,2]);
+ assert( V, [1 0 1; 0 1 1], 1e-5 );
+***** test
+ C = 2;
+ K = 3;
+ P = zeros(C,K,C,K);
+ P(1,1,1,2) = 1;
+ P(1,2,1,3) = 1;
+ P(1,3,2,2) = 1;
+ P(2,2,1,1) = 1;
+ [V ch] = qncmvisits(P);
+ assert( ch, [1 1] );
+***** test
+ C = 2;
+ K = 3;
+ P = zeros(C,K,C,K);
+ P(1,1,1,2) = 1;
+ P(1,2,1,3) = 1;
+ P(1,3,2,2) = 1;
+ P(2,2,2,1) = 1;
+ P(2,1,1,2) = 1;
+ [V ch] = qncmvisits(P);
+ assert( ch, [1 1] );
+9 tests, 9 passed, 0 known failure, 0 skipped
 [inst/ctmcisir.m]
 >>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/ctmcisir.m
 ***** test
@@ -5173,12 +5557,55 @@
  assert( max(s), 1 );
  assert( min(s), 1 );
 3 tests, 3 passed, 0 known failure, 0 skipped
+[inst/qncmnpop.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncmnpop.m
+***** test
+ H = qncmnpop( [1 2 2] );
+ assert( H, [1 1 0 0 0 0; 1 2 2 1 0 0; 1 3 5 5 3 1] );
+1 test, 1 passed, 0 known failure, 0 skipped
 [inst/qnmknode.m]
 >>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnmknode.m
 ***** test
  fail( "qnmknode( 'pippo', 1 )", "must be one" );
  fail( "qnmknode( '-/g/1-ps', 1, 1, 1)", "Invalid call" );
 1 test, 1 passed, 0 known failure, 0 skipped
+[inst/qnopen.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnopen.m
+***** test
+ # Example 34.1 p. 572
+ lambda = 3;
+ V = [16 7 8];
+ S = [0.01 0.02 0.03];
+ [U R Q X] = qnopen( lambda, S, V );
+ assert( R, [0.0192 0.0345 0.107], 1e-2 );
+ assert( U, [0.48 0.42 0.72], 1e-2 );
+***** test
+ V = [1 1; 1 1];
+ S = [1 3; 2 4];
+ lambda = [3/19 2/19];
+ [U R Q] = qnopen(lambda, S, diag( lambda / sum(lambda) ) * V);
+ assert( U(1,1), 3/19, 1e-6 );
+ assert( U(2,1), 4/19, 1e-6 );
+ assert( R(1,1), 19/12, 1e-6 );
+ assert( R(1,2), 57/2, 1e-6 );
+ assert( Q(1,1), .25, 1e-6 );
+2 tests, 2 passed, 0 known failure, 0 skipped
+[inst/ctmcfpt.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/ctmcfpt.m
+***** demo
+ Q = [ -1.0  0.9  0.1; ...
+        0.1 -1.0  0.9; ...
+        0.9  0.1 -1.0 ];
+ M = ctmcfpt(Q)
+ m = ctmcfpt(Q,1,3)
+***** test
+ N = 10;
+ Q = reshape(1:N^2,N,N);
+ Q(1:N+1:end) = 0;
+ Q -= diag(sum(Q,2));
+ M = ctmcfpt(Q);
+ assert( all(diag(M) < 10*eps) );
+1 test, 1 passed, 0 known failure, 0 skipped
 [inst/qnos.m]
 >>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnos.m
 ***** test
@@ -5205,7 +5632,7 @@
     qnos at line 136 column 30
     __test__ at line 17 column 2
     test at line 682 column 11
-    /tmp/tmp.cnygCL2JrV at line 246 column 31
+    /tmp/tmp.qG5gOSRbSf at line 358 column 31
 
 ***** test
  # Example 34.1 p. 572 Bolch et al.
@@ -5266,249 +5693,6 @@
  assert( R, [0.2 0.043 0.15 0.0625], 1e-3 );
  assert( Q, [4, 0.429 1.5 0.25], 1e-3 );
 6 tests, 6 passed, 0 known failure, 0 skipped
-[inst/erlangb.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/erlangb.m
-***** test
- fail("erlangb(1, -1)", "positive");
- fail("erlangb(-1, 1)", "positive");
- fail("erlangb(1, 0)", "positive");
- fail("erlangb(0, 1)", "positive");
- fail("erlangb('foo',1)", "positive");
- fail("erlangb(1,'bar')", "positive");
- fail("erlangb([1 1],[1 1 1])","common size");
-1 test, 1 passed, 0 known failure, 0 skipped
-[inst/dtmcexps.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/dtmcexps.m
-***** test
- P = dtmcbd([1 1 1 1], [0 0 0 0]);
- L = dtmcexps(P,[1 0 0 0 0]);
- t = dtmcmtta(P,[1 0 0 0 0]);
- assert( L, [1 1 1 1 0] );
- assert( sum(L), t );
-***** test
- P = dtmcbd(linspace(0.1,0.4,5),linspace(0.4,0.1,5));
- p0 = [1 0 0 0 0 0];
- L = dtmcexps(P,0,p0);
- assert( L, p0 );
-2 tests, 2 passed, 0 known failure, 0 skipped
-[inst/qsmmmk.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qsmmmk.m
-***** test
- lambda = mu = m = 1;
- k = 10;
- [U R Q X p0] = qsmmmk(lambda,mu,m,k);
- assert( Q, k/2, 1e-7 );
- assert( U, 1-p0, 1e-7 );
-***** test
- lambda = [1 0.8 2 9.2 0.01];
- mu = lambda + 0.17;
- k = 12;
- [U1 R1 Q1 X1] = qsmm1k(lambda,mu,k);
- [U2 R2 Q2 X2] = qsmmmk(lambda,mu,1,k);
- assert( U1, U2, 1e-5 );
- assert( R1, R2, 1e-5 );
- assert( Q1, Q2, 1e-5 );
- assert( X1, X2, 1e-5 );
- #assert( [U1 R1 Q1 X1], [U2 R2 Q2 X2], 1e-5 );
-***** test
- lambda = 0.9;
- mu = 0.75;
- k = 10;
- [U1 R1 Q1 X1 p01] = qsmmmk(lambda,mu,1,k);
- [U2 R2 Q2 X2 p02] = qsmm1k(lambda,mu,k);
- assert( [U1 R1 Q1 X1 p01], [U2 R2 Q2 X2 p02], 1e-5 );
-***** test
- lambda = 0.8;
- mu = 0.85;
- m = 3;
- k = 5;
- [U1 R1 Q1 X1 p0] = qsmmmk( lambda, mu, m, k );
- birth = lambda*ones(1,k);
- death = [ mu*linspace(1,m,m) mu*m*ones(1,k-m) ];
- q = ctmc(ctmcbd( birth, death ));
- U2 = dot( q, min( 0:k, m )/m );
- assert( U1, U2, 1e-4 );
- Q2 = dot( [0:k], q );
- assert( Q1, Q2, 1e-4 );
- assert( p0, q(1), 1e-4 );
-***** test
- # This test comes from an example I found on the web 
- lambda = 40;
- mu = 30;
- m = 3;
- k = 7;
- [U R Q X p0] = qsmmmk( lambda, mu, m, k );
- assert( p0, 0.255037, 1e-6 );
- assert( R, 0.036517, 1e-6 );
-***** test
- # This test comes from an example I found on the web 
- lambda = 50;
- mu = 10;
- m = 4;
- k = 6;
- [U R Q X p0 pk] = qsmmmk( lambda, mu, m, k );
- assert( pk, 0.293543, 1e-6 );
-***** test
- # This test comes from an example I found on the web 
- lambda = 3;
- mu = 2;
- m = 2;
- k = 5;
- [U R Q X p0 pk] = qsmmmk( lambda, mu, m, k );
- assert( p0, 0.179334, 1e-6 );
- assert( pk, 0.085113, 1e-6 );
- assert( Q, 2.00595, 1e-5 );
- assert( R-1/mu, 0.230857, 1e-6 ); # waiting time in the queue
-***** demo
- ## Given a M/M/m/K queue, compute the steady-state probability pn
- ## of having n jobs in the systen.
- lambda = 0.2;
- mu = 0.25;
- m = 5;
- K = 20;
- n = 0:10;
- pn = qsmmmk(lambda, mu, m, K, n);
- plot(n, pn, "-o", "linewidth", 2);
- xlabel("N. of jobs (n)");
- ylabel("P_n");
- title(sprintf("M/M/%d/%d system, \\lambda = %g, \\mu = %g", m, K, lambda, mu));
-7 tests, 7 passed, 0 known failure, 0 skipped
-[inst/qnomvisits.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnomvisits.m
-***** test
- fail( "qnomvisits( zeros(3,3,3), [1 1 1] )", "matrix");
-***** test
- C = 2; K = 4;
- P = zeros(C,K,C,K);
- # class 1 routing
- P(1,1,1,1) = .05;
- P(1,1,1,2) = .45;
- P(1,1,1,3) = .5;
- P(1,2,1,1) = 0.1;
- P(1,3,1,1) = 0.2;
- # class 2 routing
- P(2,1,2,1) = .01;
- P(2,1,2,3) = .5;
- P(2,1,2,4) = .49;
- P(2,3,2,1) = 0.2;
- P(2,4,2,1) = 0.16;
- lambda = [0.1 0 0 0.1 ; 0 0 0.2 0.1];
- lambda_sum = sum(lambda(:));
- V = qnomvisits(P, lambda);
- assert( all(V(:)>=0) );
- for i=1:K
-   for c=1:C
-     assert(V(c,i), lambda(c,i) / lambda_sum + sum(sum(V .* P(:,:,c,i))), 1e-5);
-   endfor
- endfor
-***** test
- # example 7.7 p. 304 Bolch et al. 
- # Note that the book uses a slightly different notation than
- # what we use here. Specifically, the book defines the routing 
- # probabilities as P(i,r,j,s) (i,j are service centers, r,s are job
- # classes) while the queueing package uses P(r,i,s,j).
- # A more serious problem arises in the definition of external arrivals.
- # The computation of V(r,i) as given in the book (called e_ir
- # in Eq 7.14) is performed in terms of P_{0, js}, defined as 
- # "the probability in an open network that a job from outside the network
- #  enters the jth node as a job of the sth class" (p. 267). This is
- # compliant with eq. 7.12 where the external class r arrival rate at center
- # i is computed as \lambda * P_{0,ir}. However, example 7.7 wrongly
- # defines P_{0,11} = P_{0,12} = 1, instead of P_{0,11} = P_{0,12} = 0.5
- # Therefore the resulting visit ratios they obtain must be divided by two.
- P = zeros(2,3,2,3);
- lambda = S = zeros(2,3);
- P(1,1,1,2) = 0.4;
- P(1,1,1,3) = 0.3;
- P(1,2,1,1) = 0.6;
- P(1,2,1,3) = 0.4;
- P(1,3,1,1) = 0.5;
- P(1,3,1,2) = 0.5;
- P(2,1,2,2) = 0.3;
- P(2,1,2,3) = 0.6;
- P(2,2,2,1) = 0.7;
- P(2,2,2,3) = 0.3;
- P(2,3,2,1) = 0.4;
- P(2,3,2,2) = 0.6;
- lambda(1,1) = lambda(2,1) = 1;
- V = qnomvisits(P,lambda);
- assert( V, [ 3.333 2.292 1.917; 10 8.049 8.415] ./ 2, 1e-3);
-3 tests, 3 passed, 0 known failure, 0 skipped
-[inst/dtmctaexps.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/dtmctaexps.m
-***** test
- P = dtmcbd([1 1 1 1], [0 0 0 0]);
- p0 = [1 0 0 0 0];
- L = dtmctaexps(P,p0);
- assert( L, [.25 .25 .25 .25 0], 10*eps );
-1 test, 1 passed, 0 known failure, 0 skipped
-[inst/qncscmva.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncscmva.m
-***** test
- N=5;
- S = [1 0.3 0.8 0.9];
- V = [1 1 1 1];
- [U1 R1 Q1 X1] = qncscmva( N, S(1:3), repmat(S(4),1,N), V );
- [U2 R2 Q2 X2] = qncsmva(N, S, V);
- assert( X1, X2, 1e-5 );
- assert( U1, U2, 1e-5 );
- assert( R1, R2, 1e-5 );
- assert( Q1, Q2, 1e-5 );
-***** test
- N=5;
- S = [1 1 1 1 1; ...
-      1 1 1 1 1; ...
-      1 1 1 1 1; ...
-      1 1/2 1/3 1/4 1/5];
- V = [1 1 1 1];
- [U1 R1 Q1 X1] = qncscmva( N, S(1:3,1), S(4,:), V );
- [U2 R2 Q2 X2] = qncsmvald(N, S, V);
- assert( U1, U2, 1e-5 );
- assert( R1, R2, 1e-5 );
- assert( Q1, Q2, 1e-5 );
- assert( X1, X2, 1e-5 );
-***** test
- N=5;
- S = [1 1 1 1 1; ...
-      1 1 1 1 1; ...
-      1 1 1 1 1; ...
-      1 1/2 1/3 1/4 1/5];
- V = [1 2 1 1];
- Z = 3;
- [U1 R1 Q1 X1] = qncscmva( N, S(1:3,1), S(4,:), V, Z );
- [U2 R2 Q2 X2] = qncsmvald(N, S, V, Z);
- assert( U1, U2, 1e-5 );
- assert( R1, R2, 1e-5 );
- assert( Q1, Q2, 1e-5 );
- assert( X1, X2, 1e-5 );
-***** demo
- maxN = 90; # Max population size
- Rmva = Rconv = Rcmva = zeros(1,maxN); # Results
- S = 4; Z = 10; m = 8;
- old = warning("query","qn:numerical-instability");
- warning("off","qn:numerical-instability");
- for N=1:maxN
-   [U R] = qncsmva(N,S,1,m,Z);		# Use MVA
-   Rmva(N) = R(1);
-   [U R] = qncsconv(N,[S Z],[1 1],[m -1]);	# Use Convolution
-   Rconv(N) = R(1);
-   if ( N > m )
-     Scmva = S ./ min(1:N,m);
-   else
-     Scmva = S ./ (1:N);
-   endif
-   [U R] = qncscmva(N,[],Scmva,1,Z);		# Use CMVA
-   Rcmva(N) = R(1);
- endfor
- warning(old.state,"qn:numerical-instability");
- plot(1:maxN, Rmva, ";MVA;", ...
-      1:maxN, Rconv, ";Convolution;", ...
-      1:maxN, Rcmva, ";CNVA;", "linewidth",2);
- xlabel("Population size (N)");
- ylabel("Response Time");
- ax=axis(); ax(3) = 0; ax(4) = 40; axis(ax);
- legend("location","northwest"); legend("boxoff");
-3 tests, 3 passed, 0 known failure, 0 skipped
 [inst/qncsvisits.m]
 >>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsvisits.m
 ***** test
@@ -5547,140 +5731,6 @@
       1 0   0  ];
  V = qncsvisits(P)
 3 tests, 3 passed, 0 known failure, 0 skipped
-[inst/erlangc.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/erlangc.m
-***** test
- fail("erlangc('foo',1)", "positive");
- fail("erlangc(1,'bar')", "positive");
-1 test, 1 passed, 0 known failure, 0 skipped
-[inst/qncsconv.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsconv.m
-***** test
- # Example 8.1 p. 318 Bolch et al.
- K=3;
- S = [ 1/0.8 1/0.6 1/0.4 ];
- m = [2 3 1];
- V = [ 1 .667 .2 ];
- [U R Q X G] = qncsconv( K, S, V, m );
- assert( G, [1 2.861 4.218 4.465], 5e-3 );
- assert( X, [0.945 0.630 0.189], 1e-3 );
- assert( U, [0.590 0.350 0.473], 1e-3 );
- assert( Q, [1.290 1.050 0.660], 1e-3 );
- assert( R, [1.366 1.667 3.496], 1e-3 );
-***** test
- # Example 8.3 p. 331 Bolch et al.
- # compare results of convolution to those of mva
- S = [ 0.02 0.2 0.4 0.6 ];
- K = 6;
- V = [ 1 0.4 0.2 0.1 ];
- [U_mva R_mva Q_mva X_mva G_mva] = qncsmva(K, S, V);
- [U_con R_con Q_con X_con G_con] = qncsconv(K, S, V);
- assert( U_mva, U_con, 1e-5 );
- assert( R_mva, R_con, 1e-5 );
- assert( Q_mva, Q_con, 1e-5 );
- assert( X_mva, X_con, 1e-5 );
- assert( G_mva, G_con, 1e-5 );
-***** test
- # Compare the results of convolution to those of mva
- S = [ 0.02 0.2 0.4 0.6 ];
- K = 6;
- V = [ 1 0.4 0.2 0.1 ];
- m = [ 1 -1 2 1 ]; # center 2 is IS
- [U_mva R_mva Q_mva X_mva] = qncsmva(K, S, V, m);
- [U_con R_con Q_con X_con G] = qncsconv(K, S, V, m );
- assert( U_mva, U_con, 1e-5 );
- assert( R_mva, R_con, 1e-5 );
- assert( Q_mva, Q_con, 1e-5 );
- assert( X_mva, X_con, 1e-5 );
-***** demo
- n = [1 2 0];
- N = sum(n); # Total population size
- S = [ 1/0.8 1/0.6 1/0.4 ];
- m = [ 2 3 1 ];
- V = [ 1 .667 .2 ];
- [U R Q X G] = qncsconv( N, S, V, m );
- p = [0 0 0]; # initialize p
- # Compute the probability to have n(k) jobs at service center k
- for k=1:3
-   p(k) = (V(k)*S(k))^n(k) / G(N+1) * ...
-          (G(N-n(k)+1) - V(k)*S(k)*G(N-n(k)) );
-   printf("Prob( n(%d) = %d )=%f\n", k, n(k), p(k) );
- endfor
-3 tests, 3 passed, 0 known failure, 0 skipped
-[inst/qnmarkov.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnmarkov.m
-***** test
- S = [5 2.5];
- P = [0 1; 1 0];
- C = [3 3];
- m = [1 1];
- [U R Q X] = qnmarkov( 3, S, C, P, m );
- assert( U, [0.9333 0.4667], 1e-4 );
- assert( X, [0.1867 0.1867], 1e-4 );
- assert( R, [12.1429 3.9286], 1e-4 );
-***** test
- S = [1/0.8 1/0.6 1/0.4];
- P = [0.6 0.3 0.1; 0.2 0.3 0.5; 0.4 0.1 0.5];
- C = [3 3 3];
- [U R Q X] = qnmarkov( 3, S, C, P );
- assert( U, [0.543 0.386 0.797], 1e-3 );
- assert( Q, [0.873 0.541 1.585], 1e-3 );
-***** xtest
- S = [2 0.9];
- C = [7 5];
- P = [0 1; 1 0];
- [U R Q X] = qnmarkov( 10, S, C, P );
- assert( Q, [6.73 3.27], 1e-3 );
-!!!!! known failure
-ASSERT errors for:  assert (Q,[6.73, 3.27],1e-3)
-
-  Location  |  Observed  |  Expected  |  Reason
-    (1)         6.4826        6.73       Abs err 0.2474 exceeds tol 0.001 by 0.2
-    (2)         3.5174        3.27       Abs err 0.2474 exceeds tol 0.001 by 0.2
-***** test
- S = [1/0.8 1/0.6 1/0.4];
- P = [(1-0.667-0.2) 0.667 0.2; 1 0 0; 1 0 0];
- m = [2 3 1];
- C = [3 3 3];
- [U R Q X] = qnmarkov( 3, S, C, P, m );
- assert( U, [0.590 0.350 0.473], 1e-3 );
- assert( Q(1:2), [1.290 1.050], 1e-3 );
-***** test
- p = 0.5; # transition prob. from S2 to S1
- mu = [1 2]; # Service rates
- C = [2 1]; # Capacities
- lambda = [0.5 0]; # arrival rate at service center 1
-
- PP = [ 0 1; p 0 ];
- [U R Q X] = qnmarkov( lambda, 1./mu, C, PP );
- ## Now we generate explicitly the infinitesimal generator matrix
- ## of the underlying MC.
- ##    00  01  10  11  20  21
- QQ = [ 0 0 lambda(1) 0 0 0; ... ## 00
-        mu(2)*(1-p) 0 mu(2)*p lambda(1) 0 0; ... ## 01
-        0 mu(1) 0 0 lambda(1) 0; ... ## 10
-        0 0 mu(2)*(1-p) 0 mu(2)*p lambda(1); ... ## 11
-        0 0 0 mu(1) 0 0; ... ## 20
-        0 0 0 0 mu(2)*(1-p) 0 ]; ## 21
- ## Complete matrix
- sum_el = sum(QQ,2);
- QQ -= diag(sum_el);
- q = ctmc(QQ);
- ## Compare results
- assert( U(1), 1-sum(q([1, 2])), 1e-5 );
- assert( U(2), 1-sum(q([1,3,5])), 1e-5 );
-***** test
- P = [0 0.5 0.5; 1 0 0; 1 0 0];
- C = [6 6 6];
- S = [1 0.8 1.8];
- N = 6;
- [U1 R1 Q1 X1] = qnclosed( N, S, qncsvisits(P) );
- [U2 R2 Q2 X2] = qnmarkov( N, S, C, P );
- assert( U1, U2, 1e-6 );
- assert( R1, R2, 1e-6 );
- assert( Q1, Q2, 1e-6 );
- assert( X1, X2, 1e-6 );
-6 tests, 5 passed, 1 known failure, 0 skipped
 [inst/qncsmvablo.m]
 >>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsmvablo.m
 ***** test
@@ -5730,101 +5780,59 @@
  assert( R1, R2, 1e-5 );
  assert( Q1, Q2, 1e-5 );
 5 tests, 5 passed, 0 known failure, 0 skipped
-[inst/qncmnpop.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncmnpop.m
-***** test
- H = qncmnpop( [1 2 2] );
- assert( H, [1 1 0 0 0 0; 1 2 2 1 0 0; 1 3 5 5 3 1] );
-1 test, 1 passed, 0 known failure, 0 skipped
-[inst/dtmc.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/dtmc.m
-***** test
- P = [0.75 0.25; 0.5 0.5];
- p = dtmc(P);
- assert( p*P, p, 1e-5 );
- assert( p, [0.6666 0.3333], 1e-4 );
-***** test
- #Example 2.11 p. 44 Bolch et al.
- P = [0.5 0.5; 0.5 0.5];
- p = dtmc(P);
- assert( p, [0.5 0.5], 1e-3 );
-***** test
- fail("dtmc( [1 1 1; 1 1 1] )", "square");
+[inst/qncsconv.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsconv.m
 ***** test
- a = 0.2;
- b = 0.8;
- P = [1-a a; b 1-b];
- plim = dtmc(P);
- p = dtmc(P, 100, [1 0]);
- assert( plim, p, 1e-5 );
+ # Example 8.1 p. 318 Bolch et al.
+ K=3;
+ S = [ 1/0.8 1/0.6 1/0.4 ];
+ m = [2 3 1];
+ V = [ 1 .667 .2 ];
+ [U R Q X G] = qncsconv( K, S, V, m );
+ assert( G, [1 2.861 4.218 4.465], 5e-3 );
+ assert( X, [0.945 0.630 0.189], 1e-3 );
+ assert( U, [0.590 0.350 0.473], 1e-3 );
+ assert( Q, [1.290 1.050 0.660], 1e-3 );
+ assert( R, [1.366 1.667 3.496], 1e-3 );
 ***** test
- P = [0 1 0 0 0; ...
-      .25 0 .75 0 0; ...
-      0 .5 0 .5 0; ...
-      0 0 .75 0 .25; ...
-      0 0 0 1 0 ];
- p = dtmc(P);
- assert( p, [.0625 .25 .375 .25 .0625], 10*eps );
+ # Example 8.3 p. 331 Bolch et al.
+ # compare results of convolution to those of mva
+ S = [ 0.02 0.2 0.4 0.6 ];
+ K = 6;
+ V = [ 1 0.4 0.2 0.1 ];
+ [U_mva R_mva Q_mva X_mva G_mva] = qncsmva(K, S, V);
+ [U_con R_con Q_con X_con G_con] = qncsconv(K, S, V);
+ assert( U_mva, U_con, 1e-5 );
+ assert( R_mva, R_con, 1e-5 );
+ assert( Q_mva, Q_con, 1e-5 );
+ assert( X_mva, X_con, 1e-5 );
+ assert( G_mva, G_con, 1e-5 );
 ***** test
- P = zeros(9,9);
- P(1,[2 4]) = 1/2;
- P(2,[1 5 3]) = 1/3;
- P(3,[2 6]) = 1/2;
- P(4,[1 5 7]) = 1/3;
- P(5,[2 4 6 8]) = 1/4;
- P(6,[3 5 9]) = 1/3;
- P(7,[4 8]) = 1/2;
- P(8,[7 5 9]) = 1/3;
- P(9,[6 8]) = 1/2;
- p = dtmc(P);
- assert( p, [1/12 1/8 1/12 1/8 1/6 1/8 1/12 1/8 1/12], 10*eps );
-***** demo
- P = zeros(9,9);
- P(1,[2 4]    ) = 1/2;
- P(2,[1 5 3]  ) = 1/3;
- P(3,[2 6]    ) = 1/2;
- P(4,[1 5 7]  ) = 1/3;
- P(5,[2 4 6 8]) = 1/4;
- P(6,[3 5 9]  ) = 1/3;
- P(7,[4 8]    ) = 1/2;
- P(8,[7 5 9]  ) = 1/3;
- P(9,[6 8]    ) = 1/2;
- p = dtmc(P);
- disp(p)
+ # Compare the results of convolution to those of mva
+ S = [ 0.02 0.2 0.4 0.6 ];
+ K = 6;
+ V = [ 1 0.4 0.2 0.1 ];
+ m = [ 1 -1 2 1 ]; # center 2 is IS
+ [U_mva R_mva Q_mva X_mva] = qncsmva(K, S, V, m);
+ [U_con R_con Q_con X_con G] = qncsconv(K, S, V, m );
+ assert( U_mva, U_con, 1e-5 );
+ assert( R_mva, R_con, 1e-5 );
+ assert( Q_mva, Q_con, 1e-5 );
+ assert( X_mva, X_con, 1e-5 );
 ***** demo
- a = 0.2;
- b = 0.15;
- P = [ 1-a a; b 1-b];
- T = 0:14;
- pp = zeros(2,length(T));
- for i=1:length(T)
-   pp(:,i) = dtmc(P,T(i),[1 0]);
+ n = [1 2 0];
+ N = sum(n); # Total population size
+ S = [ 1/0.8 1/0.6 1/0.4 ];
+ m = [ 2 3 1 ];
+ V = [ 1 .667 .2 ];
+ [U R Q X G] = qncsconv( N, S, V, m );
+ p = [0 0 0]; # initialize p
+ # Compute the probability to have n(k) jobs at service center k
+ for k=1:3
+   p(k) = (V(k)*S(k))^n(k) / G(N+1) * ...
+          (G(N-n(k)+1) - V(k)*S(k)*G(N-n(k)) );
+   printf("Prob( n(%d) = %d )=%f\n", k, n(k), p(k) );
  endfor
- ss = dtmc(P); # compute steady state probabilities
- plot( T, pp(1,:), "b+;p_0(t);", "linewidth", 2, ...
-       T, ss(1)*ones(size(T)), "b;Steady State;", ...
-       T, pp(2,:), "r+;p_1(t);", "linewidth", 2, ...
-       T, ss(2)*ones(size(T)), "r;Steady State;" );
- xlabel("Time Step");
- legend("boxoff");
-6 tests, 6 passed, 0 known failure, 0 skipped
-[inst/dtmcisir.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/dtmcisir.m
-***** test
- P = [0 .5 0; 0 0 0];
- fail( "dtmcisir(P)" );
-***** test
- P = [0 1 0; 0 .5 .5; 0 1 0];
- [r s] = dtmcisir(P);
- assert( r == 0 );
- assert( max(s), 2 );
- assert( min(s), 1 );
-***** test
- P = [.5 .5 0; .2 .3 .5; 0 .2 .8];
- [r s] = dtmcisir(P);
- assert( r == 1 );
- assert( max(s), 1 );
- assert( min(s), 1 );
 3 tests, 3 passed, 0 known failure, 0 skipped
 [inst/qncsconvld.m]
 >>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsconvld.m
@@ -5884,101 +5892,278 @@
  assert( Q, [1.290 1.050 0.660], 1e-3 );
  assert( R, [1.366 1.667 3.496], 1e-3 );
 5 tests, 5 passed, 0 known failure, 0 skipped
-[inst/ctmcexps.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/ctmcexps.m
+[inst/qnom.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnom.m
 ***** test
- Q = [-1 1; 1 -1];
- L = ctmcexps(Q,10,[1 0]);
- L = ctmcexps(Q,linspace(0,10,100),[1 0]);
+ fail( "qnom([1 1], [.9; 1.0])", "exceeded at center 1");
+ fail( "qnom([1 1], [0.9 .9; 0.9 1.0])", "exceeded at center 2");
+ #qnom([1 1], [.9; 1.0],[],2); # should not fail, M/M/2-FCFS
+ #qnom([1 1], [.9; 1.0],[],-1); # should not fail, -/G/1-PS
+ fail( "qnom(1./[2 3], [1.9 1.9 0.9; 2.9 3.0 2.9])", "exceeded at center 2");
+ #qnom(1./[2 3], [1 1.9 0.9; 0.3 3.0 1.5],[],[1 2 1]); # should not fail
 ***** test
- Q = ctmcbd( [1 2 3], [3 2 1] );
- p0 = [1 0 0 0];
- t = linspace(0,10,10);
- L1 = L2 = zeros(length(t),4);
- # compute L using the differential equation formulation
- ff = @(x,t) (x(:)'*Q+p0);
- fj = @(x,t) (Q);
- L1 = lsode( {ff, fj}, zeros(size(p0)), t );
- # compute L using ctmcexps (integral formulation)
- for i=1:length(t)
-   L2(i,:) = ctmcexps(Q,t(i),p0);
- endfor
- assert( L1, L2, 1e-5);
+ V = [1 1; 1 1];
+ S = [1 3; 2 4];
+ lambda = [3/19 2/19];
+ [U R Q X] = qnom(lambda, S, diag( lambda / sum(lambda) ) * V );
+ assert( U(1,1), 3/19, 1e-6 );
+ assert( U(2,1), 4/19, 1e-6 );
+ assert( R(1,1), 19/12, 1e-6 );
+ assert( R(1,2), 57/2, 1e-6 );
+ assert( Q(1,1), .25, 1e-6 );
+ assert( Q, R.*X, 1e-5 ); # Little's Law
+***** test
+ # example p. 138 Zahorjan et al.
+ V = [ 10 9; 5 4];
+ S = [ 1/10 1/3; 2/5 1];
+ lambda = [3/19 2/19];
+ [U R Q X] = qnom(lambda, S, diag( lambda / sum(lambda) ) * V );
+ assert( X(1,1), 1.58, 1e-2 );
+ assert( U(1,1), .158, 1e-3 );
+ assert( R(1,1), .158, 1e-3 ); # modified from the original example, as the reference above considers R as the residence time, not the response time
+ assert( Q(1,1), .25, 1e-2 );
+ assert( Q, R.*X, 1e-5 ); # Little's Law
+***** test
+ # example 7.7 p. 304 Bolch et al. Please note that the book uses the 
+ # notation P(i,r,j,s) (i,j are service centers, r,s are job
+ # classes) while the queueing package uses P(r,i,s,j)
+ P = zeros(2,3,2,3);
+ lambda = S = zeros(2,3);
+ P(1,1,1,2) = 0.4;
+ P(1,1,1,3) = 0.3;
+ P(1,2,1,1) = 0.6;
+ P(1,2,1,3) = 0.4;
+ P(1,3,1,1) = 0.5;
+ P(1,3,1,2) = 0.5;
+ P(2,1,2,2) = 0.3;
+ P(2,1,2,3) = 0.6;
+ P(2,2,2,1) = 0.7;
+ P(2,2,2,3) = 0.3;
+ P(2,3,2,1) = 0.4;
+ P(2,3,2,2) = 0.6;
+ S(1,1) = 1/8;
+ S(1,2) = 1/12;
+ S(1,3) = 1/16;
+ S(2,1) = 1/24;
+ S(2,2) = 1/32;
+ S(2,3) = 1/36;
+ lambda(1,1) = lambda(2,1) = 1;
+ V = qnomvisits(P,lambda);
+ assert( V, [ 3.333 2.292 1.917; 10 8.049 8.415] ./ 2, 1e-3);
+ [U R Q X] = qnom(sum(lambda,2), S, V);
+ assert( sum(U,1), [0.833 0.442 0.354], 1e-3 );
+ # Note: the value of K_22 (corresponding to Q(2,2)) reported in the book
+ # is 0.5. However, hand computation using the exact same formulas
+ # from the book produces a different value, 0.451
+ assert( Q, [2.5 0.342 0.186; 2.5 0.451 0.362], 1e-3 );
+***** test
+ P = zeros(2,2,2,2);
+ P(1,1,1,2) = 0.8; P(1,2,1,1) = 1;
+ P(2,1,2,2) = 0.9; P(2,2,2,1) = 1;
+ S = zeros(2,2);
+ S(1,1) = 1.5; S(1,2) = 1.2;
+ S(2,1) = 0.8; S(2,2) = 2.5;
+ lambda = zeros(2,2);
+ lambda(1,1) = 1/20;
+ lambda(2,1) = 1/30;
+ [U1 R1 Q1 X1] = qnom(lambda, S, P); # qnom_cs
+ [U2 R2 Q2 X2] = qnom(sum(lambda,2), S, qnomvisits(P,lambda)); # qnom_nocs
+ assert( U1, U2, 1e-5 );
+ assert( R1, R2, 1e-5 );
+ assert( Q1, Q2, 1e-5 );
+ assert( X1, X2, 1e-5 );
 ***** demo
- lambda = 0.5;
- N = 4;
- b = lambda*[1:N-1];
- d = zeros(size(b));
- Q = ctmcbd(b,d);
- t = linspace(0,10,100);
- p0 = zeros(1,N); p0(1)=1;
- L = zeros(length(t),N);
- for i=1:length(t)
-   L(i,:) = ctmcexps(Q,t(i),p0);
+ P = zeros(2,2,2,2);
+ lambda = zeros(2,2);
+ S = zeros(2,2);
+ P(1,1,2,1) = P(1,2,2,1) = 0.2;
+ P(1,1,2,2) = P(2,2,2,2) = 0.8;
+ S(1,1) = S(1,2) = 0.1;
+ S(2,1) = S(2,2) = 0.05;
+ rr = 1:100;
+ Xk = zeros(2,length(rr));
+ for r=rr
+   lambda(1,1) = lambda(1,2) = 1/r;
+   [U R Q X] = qnom(lambda,S,P);
+   Xk(:,r) = sum(X,1)';
  endfor
- plot( t, L(:,1), ";State 1;", "linewidth", 2, ...
-       t, L(:,2), ";State 2;", "linewidth", 2, ...
-       t, L(:,3), ";State 3;", "linewidth", 2, ...
-       t, L(:,4), ";State 4;", "linewidth", 2 );
- legend("location","northwest"); legend("boxoff");
- xlabel("Time");
- ylabel("Expected sojourn time");
-***** demo
- lambda = 0.5;
- N = 4;
- b = lambda*[1:N-1];
- d = zeros(size(b));
- Q = ctmcbd(b,d);
- p0 = zeros(1,N); p0(1)=1;
- L = ctmcexps(Q,p0);
- disp(L);
-2 tests, 2 passed, 0 known failure, 0 skipped
-[inst/qncsaba.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsaba.m
+ plot(rr,Xk(1,:),";Server 1;","linewidth",2, ...
+      rr,Xk(2,:),";Server 2;","linewidth",2);
+ legend("boxoff");
+ xlabel("Class 1 interarrival time");
+ ylabel("Throughput");
+5 tests, 5 passed, 0 known failure, 0 skipped
+[inst/qncmmva.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncmmva.m
 ***** test
- fail("qncsaba(-1,0)", "N must be");
- fail("qncsaba(1,[])", "nonempty");
- fail("qncsaba(1,[-1 2])", "nonnegative");
- fail("qncsaba(1,[1 2],[1 2 3])", "incompatible size");
- fail("qncsaba(1,[1 2 3],[1 2 -1])", "nonnegative");
- fail("qncsaba(1,[1 2 3],[1 2 3],[1 2])", "incompatible size");
- fail("qncsaba(1,[1 2 3],[1 2 3],[1 2 1])", "not supported");
- fail("qncsaba(1,[1 2 3],[1 2 3],[1 1 1],-1)", "nonnegative");
- fail("qncsaba(1,[1 2 3],[1 2 3],1,[0 0])", "scalar");
+ S = [1 1 2; 1 1 1];
+ V = [1 1 1; 1 1 1];
+ N = [1 1];
+ m = [1 1 2];
+ fail( "qncmmva(N)" );
+ fail( "qncmmva(N,S,V,m)", "independent" );
+ S = [0 0 0; 1 1 1];
+ fail( "qncmmva(N,S,V,m)", "must contain at least" );
+ S = [1 2 3; 1 2 3];
+ N = [1 1];
+ V = zeros(3,2,3);
+ fail( "qncmmva(N,S,V)", "size mismatch" );
+ fail( "qncmmva([0.3 1], [1 2; 3 4])", "integer");
+ fail( "qncmmva([-1 0], [1 2; 3 4])", "nonnegative");
+warning: the '.+' operator was deprecated in version 7 and will not be allowed in a future version of Octave; please use '+' instead; near line 430 of file '/build/reproducible-path/octave-queueing-1.2.7/debian/octave-queueing/usr/share/octave/packages/queueing-1.2.7/qncmmva.m'
+warning: called from
+    print_usage at line 62 column 18
+    qncmmva at line 204 column 5
+    fail at line 130 column 7
+    __test__ at line 7 column 2
+    test at line 682 column 11
+    /tmp/tmp.qG5gOSRbSf at line 406 column 31
+
 ***** test
- N = 20;
- S = [ 4.6*2 8 ];
- Z = 120;
- [X_l X_u R_l R_u] = qncsaba(N, S, ones(size(S)), ones(size(S)), Z);
- assert( [X_u R_l], [0.109 64], 1e-3 );
-2 tests, 2 passed, 0 known failure, 0 skipped
-[inst/qncmvisits.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncmvisits.m
+ S = [1 1 1; 1 1 1];
+ N = [0 0];
+ [U R Q X] = qncmmva(N, S);
+ assert( U, 0*S );
+ assert( R, 0*S );
+ assert( Q, 0*S );
+ assert( X, 0*S );
 ***** test
-
- ## Closed, multiclass network
-
- C = 2; K = 3; 
+ S = [1 1 1; 1 1 1];
+ V = [1 1 1; 1 1 1];
+ N = [0 0];
+ m = [2 2 2];
+ [U R Q X] = qncmmva(N, S, V, m);
+ assert( U, 0*S );
+ assert( R, 0*S );
+ assert( Q, 0*S );
+ assert( X, 0*S );
+***** test
+ S = [ 1/10 1/3; 2/5 1 ];
+ V = [ 10 9; 5 4 ];
+ N = [ 1 1 ];
+ [U R Q X] = qncmmva(N,S,V);
+ assert( Q, [ 4/19 15/19; 5/19 14/19 ], 1e-3 );
+ assert( R .* V, [ 4/3 5; 5/2 7 ], 1e-3 );
+ assert( diag( X ./ V )', [ 3/19 2/19 ], 1e-3 );
+ assert( all(U(:)<=1) );
+ assert( Q, R.*X, 1e-5 ); # Little's Law
+***** test
+ S = [0.02 0.2 0.4 0.6];
+ V = [1 0.4 0.2 0.1];
+ N = [6];
+ [U R Q X] = qncmmva( N, S, V );
+ assert( Q, [0.244 2.261 2.261 1.234], 1e-3 );
+ assert( R, [0.025 0.570 1.140 1.244], 1e-3 );
+ assert( X, [9.920 3.968 1.984 0.992], 1e-3 );
+ assert( U, [0.198 0.794 0.794 0.595], 1e-3 );
+***** test
+ S = [1 0 .025; 0 15 .5];
+ V = [1 0 1; 0 1 1];
+ N = [2 1];
+ m = [-1 -1 1];
+ [U R Q X] = qncmmva(N,S,V,m);
+ assert( R(1,1), 1, 1e-3 );
+ assert( R(2,2), 15, 1e-3 );
+ assert( R(1,3), .027, 1e-3 );
+ assert( R(2,3), .525, 1e-3  );
+ assert( X(1,1)+X(1,2), 1.949, 1e-3 );
+ assert( X(2,1)+X(2,2), 0.064, 1e-3 );
+ assert( sum(Q,1), [1.949, .966, .085], 1e-3 );
+ assert( all(U(:,3)<=1) );
+ assert( Q, R.*X, 1e-5 ); # Little's Law
+***** test
+ S = [1 0 .025; 0 15 .5];
+ V = [1 0 1; 0 1 1];
+ N = [15 5];
+ m = [-1 -1 1];
+ [U R Q X] = qncmmva(N,S,V,m);
+ # FIXME: I replaced 14.3->14.323
+ assert( U, [14.323 0 .358; 0 4.707 .157], 1e-3 );
+ # FIXME: I replaced 14.3->14.323
+ assert( X, [14.323 0 14.323; 0 .314 .314 ], 1e-3 );
+ # FIXME: I replaced 14.3->14.323
+ assert( Q, [14.323 0 .677; 0 4.707 .293 ], 1e-3 );
+ assert( R, [1 0 .047; 0 15 .934 ], 1e-3 );
+ assert( Q, R.*X, 1e-5 ); # Little's Law
+***** test
+ S = [ 0.2 0.4 1; 0.2 0.6 2 ];
+ V = [ 1 0.6 0.4; 1 0.3 0.7 ];
+ N = [ 2 1 ];
+ m = [ 2 1 -1 ];
+ [U R Q X] = qncmmva(N,S,V,m);
+ assert( Q, [ 0.428 0.726 0.845; 0.108 0.158 0.734 ], 1e-3 );
+ assert( X(1,1), 2.113, 1e-3 ); # CHECK
+ assert( X(2,1), 0.524, 1e-3 ); # CHECK
+ assert( all(U(:)<=1) );
+ assert( Q, R.*X, 1e-5 ); # Little's Law
+***** test
+ C = 2; # two classes
+ K = 4; # four servers
+ S = V = zeros(C,K);
+ S(1,:) = linspace(1,2,K);
+ S(2,:) = linspace(2,3,K);
+ V(1,:) = linspace(4,1,K);
+ V(2,:) = linspace(6,3,K);
+ N = [10 0]; # class 2 has no customers
+ [U1 R1 Q1 X1] = qncmmva(N,S,V);
+ [U2 R2 Q2 X2] = qncsmva(N(1),S(1,:),V(1,:));
+ assert( U1(1,:), U2, 1e-5 );
+ assert( R1(1,:), R2, 1e-5 );
+ assert( Q1(1,:), Q2, 1e-5 );
+ assert( X1(1,:), X2, 1e-5 );
+***** test
+ Z = [1 15];
+ V = [1; 1];
+ S = [.025; .5];
+ N = [15; 5];
+ [U R Q X] = qncmmva(N, S, V, 1, Z);
+ assert( U, [.358; .157], 1e-3 );
+ assert( Q, [.677; .293], 1e-3 );
+ assert( X, [14.323; .314], 1e-3 ); ## NOTE: X(1,1) = 14.3 in Schwetman
+ assert( R, [.047; .934], 1e-3 );
+ assert( Q, R.*X, 1e-5 ); # Little's Law
+***** test
+ C = 2;
+ K = 6;
  P = zeros(C,K,C,K);
- P(1,1,1,2) = 1;
- P(1,2,1,1) = 1;
- P(2,1,2,3) = 1;
- P(2,3,2,1) = 1;
+ P(1,1,1,2) = P(2,1,2,2) = 1;
+ P(1,2,1,3) = P(1,2,1,4) = P(1,2,1,5) = P(1,2,1,6) = .25;
+ P(2,2,2,3) = P(2,2,2,4) = P(2,2,2,5) = P(2,2,2,6) = .25;
+ P(1,3,1,1) = P(1,4,1,1) = P(1,5,1,1) = P(1,6,1,1) = .9;
+ P(1,3,1,2) = P(1,4,1,2) = P(1,5,1,2) = P(1,6,1,2) = .1;
+ P(2,3,2,1) = P(2,4,2,1) = P(2,5,2,1) = P(2,6,2,1) = .05;
+ P(2,3,2,2) = P(2,4,2,2) = P(2,5,2,2) = P(2,6,2,2) = .95;
+ N = [40 4];
+ S = [ 5.0 .010 .035 .035 .035 .035; ...
+      10.0 .100 .035 .035 .035 .035 ];
  V = qncmvisits(P);
- for c=1:C
-   for k=1:K
-     assert(V(c,k), sum(sum(V .* P(:,:,c,k))), 1e-5);
-   endfor
- endfor
+ [U R Q X] = qncmmva(N, S, V, [-1 1 1 1 1 1]);
+ # FIXME: The results below were computed with JMVA; the numbers
+ # in the paper are different (wrong?!?)!!
+ assert( U, [39.457941 0.087684 0.076724 0.076724 0.076724 0.076724; ...
+              2.772704 0.554541 0.048522 0.048522 0.048522 0.048522 ], 1e-5 );
+ assert( R.*V, [5 0.024363 0.011081 0.011081 0.011081 0.011081; ...
+                10 3.636155 0.197549 0.197549 0.197549 0.197549 ], 1e-5 );
+ assert( Q(:,1), [39.457941 2.772704]', 1e-5 );
+ assert( Q(:,2), [0.192262 1.008198]', 1e-5 );
+ assert( Q(:,3), [0.087449 0.054775]', 1e-5 );
+ assert( Q(:,4), Q(:,5), 1e-5 );
+ assert( Q(:,5), Q(:,6), 1e-5 );
+ assert( X(:,1), [7.891588 0.277270]', 1e-5 );
+ assert( X(:,2), [8.768431 5.545407]', 1e-5 );
+ assert( X(:,3), [2.192108 1.386352]', 1e-5 );
+ assert( X(:,4), X(:,5), 1e-5 );
+ assert( X(:,5), X(:,6), 1e-5 );
+ assert( Q, R.*X, 1e-5 ); # Little's Law
 ***** test
-
- ## Test multiclass network. Example from Schwetman (figure 7, page 9 of
- ## http://docs.lib.purdue.edu/cstech/259/
- ## "Testing network-of-queues software, technical report CSD-TR 330,
- ## Purdue University).
-
+ C = 2; # two classes
+ K = 4; # four servers
  C = 2; K = 4;
  P = zeros(C,K,C,K);
+ S = zeros(C,K);
+
+ # Routing
+
  # class 1 routing
  P(1,1,1,1) = .05;
  P(1,1,1,2) = .45;
@@ -5991,24 +6176,53 @@
  P(2,1,2,4) = .49;
  P(2,3,2,1) = 1;
  P(2,4,2,1) = 1;
- V = qncmvisits(P);
- for c=1:C
-   for i=1:K
-     assert(V(c,i), sum(sum(V .* P(:,:,c,i))), 1e-5);
-   endfor
- endfor
-***** test
 
- ## Network with class switching.
- ## This is the example in figure 9 of
- ## Schwetman, "Implementing the Mean Value Analysis
- ## Algorithm fort the solution of Queueing Network Models", Technical
- ## Report CSD-TR-355, Department of Computer Science, Purdue Univrsity,
- ## Feb 15, 1982, http://docs.lib.purdue.edu/cstech/286/
+ # Compute visits
+
+ V = qncmvisits(P);
 
+ # Define population and service times
+ 
+ N = [3 2];
+ S = [0.01 0.09 0.10 0.08; ...
+      0.05 0.09 0.10 0.08];
+ [U1 R1 Q1 X1] = qncmmva(N,S,V); # this invokes __qncmmva_nocs
+ [U2 R2 Q2 X2] = qncmmva(N,S,P); # this invokes __qncmmva_cs
+ assert( U2, U1, 1e-5 );
+ assert( R2, R1, 1e-5 );
+ assert( Q2, Q1, 1e-5 );
+ assert( X2, X1, 1e-5 );
+***** test
+ S = [1 0 .025; 0 15 .5];
+ V = [1 0 1; 0 1 1];
+ N = [15 5];
+ m = [-1 -1 1];
+ [U R Q X] = qncmmva(N,S,V,m);
+ assert( U, [14.323 0 .358; 0 4.707 .157], 1e-3 );
+ assert( R, [1.0 0 .047; 0 15 .934], 1e-3 );
+ assert( Q, [14.323 0 .677; 0 4.707 .293], 1e-3 );
+ assert( X, [14.323 0 14.323; 0 .314 .314], 1e-3 );
+***** test
+ S = [.025 1 15; .5 1 15 ];
+ P = zeros(2,3,2,3);
+ P(1,1,1,2) = P(1,2,1,1) = 1;
+ P(2,1,2,3) = P(2,3,2,1) = 1;
+ N = [15 5];
+ m = [1 -1 -1];
+ r = [1 1]; # reference station is station 1
+ [U R Q X] = qncmmva(N,S,P,r,m);
+ # FIXME: I replaced 14.3->14.323
+ assert( U, [0.358 14.323 0; 0.156 0 4.707], 1e-3 );
+ # FIXME: I replaced 14.3->14.323
+ assert( X, [14.323 14.3230 0; .314 0 .314 ], 1e-3 );
+ # FIXME: I replaced 14.3->14.323
+ assert( Q, [.677 14.323 0; .293 0 4.707], 1e-3 );
+ assert( R, [.047 1 15.0; .934 1 15.0], 1e-3 );
+ assert( Q, R.*X, 1e-5 ); # Little's Law
+***** test
  C = 2; K = 3;
  S = [.01 .07 .10; ...
-      .05 0.7 .10 ];
+      .05 .07 .10 ];
  P = zeros(C,K,C,K);
  P(1,1,1,2) = .7;
  P(1,1,1,3) = .2;
@@ -6016,380 +6230,135 @@
  P(2,1,2,2) = .3;
  P(2,1,2,3) = .5;
  P(2,1,1,1) = .2;
- P(1,2,1,1) = P(1,3,1,1) = 1;
- P(2,2,2,1) = P(2,3,2,1) = 1;
+ P(1,2,1,1) = P(2,2,2,1) = 1;
+ P(1,3,1,1) = P(2,3,2,1) = 1;
  N = [3 0];
- V = qncmvisits(P);
- VV = [10 7 2; 5 1.5 2.5]; # result given in Schwetman; our function computes something different, but that's ok since visit counts are actually ratios
- assert( V ./ repmat(V(:,1),1,K), VV ./ repmat(VV(:,1),1,K), 1e-5 );
-***** test
-
- ## two disjoint classes: must produce two disjoing chains
-
- C = 2; K = 3;
- P = zeros(C,K,C,K);
- P(1,1,1,2) = 1;
- P(1,2,1,1) = 1;
- P(2,1,2,3) = 1;
- P(2,3,2,1) = 1;
- [nc r] = qncmvisits(P);
- assert( r(1) != r(2) );
-***** test
-
- ## two classes, one chain
-
- C = 2; K = 3;
- P = zeros(C,K,C,K);
- P(1,1,1,2) = .5;
- P(1,2,2,1) = 1;
- P(2,1,2,3) = .5;
- P(2,3,1,1) = 1;
- [nc r] = qncmvisits(P);
- assert( r(1) == r(2) );
+ [U R Q X] = qncmmva(N, S, P);
+ assert( R, [.015 .133 .163; .073 .133 .163], 1e-3 );
+ assert( X, [12.609 8.826 2.522; 6.304 1.891 3.152], 1e-3 );
+ assert( Q, [.185 1.175 .412; .462 .252 .515], 1e-3 );
+ assert( U, [.126 .618 .252; .315 .132 .315], 1e-3 );
 ***** test
- 
- ## a "Moebius strip". Note that this configuration is invalid, and
- ## therefore our algorithm must raise an error. This is because this
- ## network has two chains, but both chains contain both classes
-
- C = 2; K = 2;
- P = zeros(C,K,C,K);
- P(1,1,2,2) = 1;
- P(2,2,1,1) = 1;
- P(2,1,1,2) = 1;
- P(1,2,2,1) = 1;
- fail( "qncmvisits(P)", "different");
+ V = [ 1.00 0.45 0.50 0.00; ...
+       1.00 0.00 0.50 0.49 ];
+ N = [3 2];
+ S = [0.01 0.09 0.10 0.08; ...
+      0.05 0.09 0.10 0.08];
+ [U R Q X] = qncmmva(N, S, V);
+ assert( U, [ 0.1215 0.4921 0.6075 0.0000; ...
+              0.3433 0.0000 0.3433 0.2691 ], 1e-4 );
+ assert( Q, [ 0.2131 0.7539 2.0328 0.0000; ...
+              0.5011 0.0000 1.1839 0.3149 ], 1e-4 );
+ assert( R.*V, [0.0175 0.0620 0.1672 0.0000; ...
+                0.0729 0.0000 0.1724 0.0458 ], 1e-4 );
+ assert( X, [12.1517 5.4682 6.0758 0.0000; ...
+              6.8669 0.0000 3.4334 3.3648 ], 1e-4 );
 ***** test
-
- ## Network with two classes representing independent chains.
- ## This is example in figure 8 of
- ## Schwetman, "Implementing the Mean Value Analysis
- ## Algorithm fort the solution of Queueing Network Models", Technical
- ## Report CSD-TR-355, Department of Computer Science, Purdue Univrsity,
- ## Feb 15, 1982, http://docs.lib.purdue.edu/cstech/286/
- 
- C = 2; K = 2;
- P = zeros(C,K,C,K);
+ N = [1   1   1];
+ S = [0.20000   0.02000   0.00000   0.00000;
+      0.20000   0.00000   0.02000   0.00000;
+      0.20000   0.00000   0.00000   0.02000];
+ V = [1   1   0   0;
+      1   0   1   0;
+      1   0   0   1];
+ [U R Q X] = qncmmva(N,S,V);
+ assert( Q(1,2), Q(2,3), 1e-5);
+ assert( Q(2,3), Q(3,4), 1e-5);
+ assert( abs(max(Q(:,1)) - min(Q(:,1))) < 1e-5 );
+***** demo
+ Ntot = 100; # total population size
+ b = linspace(0.1,0.9,10); # fractions of class-1 requests
+ S = [20 80 31 14 23 12; ...
+      90 30 33 20 14 7];
+ V = ones(size(S));
+ X1 = X1 = XX = zeros(size(b));
+ R1 = R2 = RR = zeros(size(b));
+ for i=1:length(b)
+   N = [fix(b(i)*Ntot) Ntot-fix(b(i)*Ntot)];
+   # printf("[%3d %3d]\n", N(1), N(2) );
+   [U R Q X] = qncmmva( N, S, V );
+   X1(i) = X(1,1) / V(1,1);
+   X2(i) = X(2,1) / V(2,1);
+   XX(i) = X1(i) + X2(i);
+   R1(i) = dot(R(1,:), V(1,:));
+   R2(i) = dot(R(2,:), V(2,:));
+   RR(i) = Ntot / XX(i);
+ endfor
+ subplot(2,1,1);
+ plot(b, X1, ";Class 1;", "linewidth", 2, ...
+      b, X2, ";Class 2;", "linewidth", 2, ...
+      b, XX, ";System;", "linewidth", 2 );
+ legend("location","south"); legend("boxoff");
+ ylabel("Throughput");
+ subplot(2,1,2);
+ plot(b, R1, ";Class 1;", "linewidth", 2, ...
+      b, R2, ";Class 2;", "linewidth", 2, ...
+      b, RR, ";System;", "linewidth", 2 );
+ legend("location","south"); legend("boxoff");
+ xlabel("Population mix \\beta for Class 1");
+ ylabel("Resp. Time");
+***** demo
+ S = [1 0 .025; 0 15 .5];
+ P = zeros(2,3,2,3);
  P(1,1,1,3) = P(1,3,1,1) = 1;
  P(2,2,2,3) = P(2,3,2,2) = 1;
- V = qncmvisits(P,[1,2]);
- assert( V, [1 0 1; 0 1 1], 1e-5 );
-***** test
- C = 2;
- K = 3;
- P = zeros(C,K,C,K);
- P(1,1,1,2) = 1;
- P(1,2,1,3) = 1;
- P(1,3,2,2) = 1;
- P(2,2,1,1) = 1;
- [V ch] = qncmvisits(P);
- assert( ch, [1 1] );
-***** test
- C = 2;
- K = 3;
+ V = qncmvisits(P,[3 3]); # reference station is station 3
+ N = [15 5];
+ m = [-1 -1 1];
+ [U R Q X] = qncmmva(N,S,V,m)
+***** demo
+ C = 2; K = 3;
+ S = [.01 .07 .10; ...
+      .05 .07 .10 ];
  P = zeros(C,K,C,K);
- P(1,1,1,2) = 1;
- P(1,2,1,3) = 1;
- P(1,3,2,2) = 1;
- P(2,2,2,1) = 1;
- P(2,1,1,2) = 1;
- [V ch] = qncmvisits(P);
- assert( ch, [1 1] );
-9 tests, 9 passed, 0 known failure, 0 skipped
-[inst/qncsgb.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsgb.m
-***** test
- fail( "qncsgb( 1, [] )", "vector" );
- fail( "qncsgb( 1, [0 -1])", "nonnegative" );
- fail( "qncsgb( 0, [1 2] )", "> 0" );
- fail( "qncsgb( -1, [1 2])", "nonnegative" );
- fail( "qncsgb( 1, [1 2],1,[1 -1])", "single server" );
-warning: the '.-' operator was deprecated in version 7 and will not be allowed in a future version of Octave; please use '-' instead; near line 182 of file '/build/reproducible-path/octave-queueing-1.2.7/debian/octave-queueing/usr/share/octave/packages/queueing-1.2.7/qncsgb.m'
-warning: called from
-    fail at line 130 column 7
-    __test__ at line 3 column 2
-    test at line 682 column 11
-    /tmp/tmp.cnygCL2JrV at line 398 column 31
-
-warning: the '.-' operator was deprecated in version 7 and will not be allowed in a future version of Octave; please use '-' instead; near line 184 of file '/build/reproducible-path/octave-queueing-1.2.7/debian/octave-queueing/usr/share/octave/packages/queueing-1.2.7/qncsgb.m'
-warning: called from
-    fail at line 130 column 7
-    __test__ at line 3 column 2
-    test at line 682 column 11
-    /tmp/tmp.cnygCL2JrV at line 398 column 31
-
-***** # shared test function
-***** function test_gb( D, expected, Z=0 )
- for i=1:rows(expected)
-   N = expected(i,1);
-   [X_lower X_upper Q_lower Q_upper] = qncsgb(N,D,1,1,Z);
-   X_exp_lower = expected(i,2);
-   X_exp_upper = expected(i,3);
-   assert( [N X_lower X_upper], [N X_exp_lower X_exp_upper], 1e-4 )
- endfor
-***** xtest
- # table IV
- D = [ 0.1 0.1 0.09 0.08 ];
- #            N  X_lower  X_upper
- expected = [ 2  4.3040   4.3174; ...
-              5  6.6859   6.7524; ...
-              10 8.1521   8.2690; ...
-              20 9.0947   9.2431; ...
-              80 9.8233   9.8765 ];
- test_gb(D, expected);
-!!!!! known failure
-ASSERT errors for:  assert ([N, X_lower, X_upper],[N, X_exp_lower, X_exp_upper],1e-4)
-
-  Location  |  Observed  |  Expected  |  Reason
-    (3)         4.3813       4.3174      Abs err 0.063891 exceeds tol 0.0001 by 0.06
-***** xtest
- # table V
- D = [ 0.1 0.1 0.09 0.08 ];
- Z = 1;
- #            N  X_lower  X_upper
- expected = [ 2  1.4319   1.5195; ...
-              5  3.3432   3.5582; ...
-              10 5.7569   6.1410; ...
-              20 8.0856   8.6467; ...
-              80 9.7147   9.8594];
- test_gb(D, expected, Z);
-!!!!! known failure
-ASSERT errors for:  assert ([N, X_lower, X_upper],[N, X_exp_lower, X_exp_upper],1e-4)
-
-  Location  |  Observed  |  Expected  |  Reason
-    (3)         1.5219       1.5195      Abs err 0.0023822 exceeds tol 0.0001 by 0.002
-***** test
- P = [0 0.3 0.7; 1 0 0; 1 0 0];
- S = [1 0.6 0.2];
- m = ones(1,3);
- V = qncsvisits(P);
- Z = 2;
- Nmax = 20;
- tol = 1e-5; # compensate for numerical errors
- ## Test case with Z>0
- for n=1:Nmax
-   [X_gb_lower X_gb_upper NC NC Q_gb_lower Q_gb_upper] = qncsgb(n, S.*V, 1, 1, Z);
-   [U R Q X] = qnclosed( n, S, V, m, Z );
-   X_mva = X(1)/V(1);
-   assert( X_gb_lower <= X_mva+tol );
-   assert( X_gb_upper >= X_mva-tol );
-   assert( Q_gb_lower <= Q+tol ); # compensate for numerical errors
-   assert( Q_gb_upper >= Q-tol ); # compensate for numerical errors
- endfor
-***** test
- P = [0 0.3 0.7; 1 0 0; 1 0 0];
- S = [1 0.6 0.2];
- V = qncsvisits(P);
- Nmax = 20;
- tol = 1e-5; # compensate for numerical errors
-
- ## Test case with Z=0
- for n=1:Nmax
-   [X_gb_lower X_gb_upper NC NC Q_gb_lower Q_gb_upper] = qncsgb(n, S.*V);
-   [U R Q X] = qnclosed( n, S, V );
-   X_mva = X(1)/V(1);
-   assert( X_gb_lower <= X_mva+tol );
-   assert( X_gb_upper >= X_mva-tol );
-   assert( Q_gb_lower <= Q+tol );
-   assert( Q_gb_upper >= Q-tol );
- endfor
-5 tests, 3 passed, 2 known failures, 0 skipped
-[inst/qncsbsb.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncsbsb.m
-***** test
- fail("qncsbsb(-1,0)", "N must be");
- fail("qncsbsb(1,[])", "nonempty");
- fail("qncsbsb(1,[-1 2])", "nonnegative");
- fail("qncsbsb(1,[1 2],[1 2 3])", "incompatible size");
- fail("qncsbsb(1,[1 2 3],[1 2 3],[1 2])", "incompatible size");
- fail("qncsbsb(1,[1 2 3],[1 2 3],[1 2 1])", "M/M/1 servers");
- fail("qncsbsb(1,[1 2 3],[1 2 3],[1 1 1],-1)", "nonnegative");
- fail("qncsbsb(1,[1 2 3],[1 2 3],[1 1 1],[0 0])", "scalar");
-1 test, 1 passed, 0 known failure, 0 skipped
-[inst/qnsolve.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnsolve.m
-***** test
- # Example 8.7 p. 349 Bolch et al.
- N = 3;
- Q1 = qnmknode( "m/m/m-fcfs", .5, 2 );
- Q2 = qnmknode( "m/m/m-fcfs", 1/1.667 );
- Q3 = qnmknode( "m/m/m-fcfs", 1/1.25 );
- Q4 = qnmknode( "m/m/m-fcfs", 1./[1 2 3] );
- V = [ 1 .5 .5 1 ];
- [U R Q X] = qnsolve("closed",N, { Q1, Q2, Q3, Q4 }, V);
- assert( Q, [0.624 0.473 0.686 1.217], 1e-3 );
- assert( R, [0.512 0.776 1.127 1], 1e-3 );
-warning: the '.+' operator was deprecated in version 7 and will not be allowed in a future version of Octave; please use '+' instead; near line 424 of file '/build/reproducible-path/octave-queueing-1.2.7/debian/octave-queueing/usr/share/octave/packages/queueing-1.2.7/qnsolve.m'
-warning: called from
-    __test__ at line 10 column 15
-    test at line 682 column 11
-    /tmp/tmp.cnygCL2JrV at line 414 column 31
-
-***** test
- # Example 8.7 p. 349 Bolch et al.
- N = 3;
- Q1 = qnmknode( "m/m/m-fcfs", 1/2, 2 );
- Q2 = qnmknode( "m/m/m-fcfs", 1/1.667 );
- Q3 = qnmknode( "m/m/m-fcfs", 1/1.25 );
- Q4 = qnmknode( "-/g/inf", 1 );
- V = [ 1 .5 .5 1 ];
- [U R Q X] = qnsolve("closed",N, { Q1, Q2, Q3, Q4 }, V);
- assert( Q, [0.624 0.473 0.686 1.217], 1e-3 );
- assert( R, [0.512 0.776 1.127 1], 1e-3 );
-***** test
- # Example 8.4 p. 333 Bolch et al.
- N = 3;
- Q1 = qnmknode( "m/m/m-fcfs", .5, 2 );
- Q2 = qnmknode( "m/m/m-fcfs", .6 );
- Q3 = qnmknode( "m/m/m-fcfs", .8 );
- Q4 = qnmknode( "-/g/inf", 1 );
- V = [ 1 .5 .5 1 ];
- [U R Q X] = qnsolve("closed",N, { Q1, Q2, Q3, Q4 }, V);
- assert( U(1:3), [.304 .365 .487], 1e-3 );
- assert( X, [1.218 0.609 0.609 1.218], 1e-3 );
-***** test
- # Same as above, with center 1 replaced with a load-dependent service center
- N = 3;
- Q1 = qnmknode( "m/m/m-fcfs", [.5 .25 .25] );
- Q2 = qnmknode( "m/m/m-fcfs", .6 );
- Q3 = qnmknode( "m/m/m-fcfs", .8 );
- Q4 = qnmknode( "m/m/m-fcfs", [1 1/2 1/3] );
- V = [ 1 .5 .5 1 ];
- [U R Q X] = qnsolve("closed",N, { Q1, Q2, Q3, Q4 }, V);
- assert( U(2:3), [.365 .487], 1e-3 ); ## NOTE: Utilization U(1) is computed differently from M/M/m nodes and load-dependent M/M/1 nodes
- assert( X, [1.218 0.609 0.609 1.218], 1e-3 );
-***** test
- # Example 7.4 p. 287 Bolch et al.
- QQ = { qnmknode( "m/m/m-fcfs", 0.04 ), ...
-        qnmknode( "m/m/m-fcfs", 0.03 ), ...
-        qnmknode( "m/m/m-fcfs", 0.06 ), ...
-        qnmknode( "m/m/m-fcfs", 0.05 ) };
- P = [ 0 0.5 0.5 0; 1 0 0 0; 0.6 0 0 0; 1 0 0 0 ];
- lambda = [0 0 0 4];
- [U R Q X] = qnsolve("open", sum(lambda), QQ, qnosvisits(P,lambda) );
- assert( X, [20 10 10 4], 1e-4 );
- assert( U, [0.8 0.3 0.6 0.2], 1e-2 );
- assert( R, [0.2 0.043 0.15 0.0625], 1e-3 );
- assert( Q, [4, 0.429 1.5 0.25], 1e-3 );
-***** test
- V = [1 1; 1 1];
- Q1 = qnmknode( "m/m/m-fcfs", [1;2] );
- Q2 = qnmknode( "m/m/m-fcfs", [3;4] );
- lambda = [3/19 2/19];
- [U R Q] = qnsolve("open", lambda, { Q1, Q2 }, diag( lambda / sum(lambda) ) * V);
- assert( U(1,1), 3/19, 1e-6 );
- assert( U(2,1), 4/19, 1e-6 );
- assert( R(1,1), 19/12, 1e-6 );
- assert( R(1,2), 57/2, 1e-6 );
- assert( Q(1,1), .25, 1e-6 );
+ P(1,1,1,2) = .7; P(1,1,1,3) = .2; P(1,1,2,1) = .1;
+ P(2,1,2,2) = .3; P(2,1,2,3) = .5; P(2,1,1,1) = .2;
+ P(1,2,1,1) = P(2,2,2,1) = 1;
+ P(1,3,1,1) = P(2,3,2,1) = 1;
+ N = [3 0];
+ [U R Q X] = qncmmva(N, S, P)
+17 tests, 17 passed, 0 known failure, 0 skipped
+[inst/dtmcchkP.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/dtmcchkP.m
 ***** test
- QQ = { qnmknode( "m/m/m-fcfs", [0.2; 0.2], 2 ), ...
-        qnmknode( "-/g/1-ps", [0.4; 0.6] ), ...
-        qnmknode( "-/g/inf", [1; 2] ) };
- V = [ 1 0.6 0.4; 1 0.3 0.7 ];
- N = [ 2 1 ];
- [U R Q X] = qnsolve( "closed", N, QQ, V );
- assert( Q, [ 0.428 0.726 0.845; 0.108 0.158 0.734 ], 1e-3 );
- assert( X(1,1), 2.113, 1e-3 ); # CHECK
- assert( X(2,1), 0.524, 1e-3 ); # CHECK
- assert( all( all(U(:,[1,2])<=1) ) );
+ [r err] = dtmcchkP( [1 1 1; 1 1 1] );
+ assert( r, 0 );
+ assert( index(err, "square") > 0 );
 ***** test
- QQ = { qnmknode( "m/m/m-fcfs", [0.2 0.1 0.1; 0.2 0.1 0.1] ), ...
-        qnmknode( "-/g/1-ps", [0.4; 0.6] ), ...
-        qnmknode( "-/g/inf", [1; 2] ) };
- V = [ 1 0.6 0.4; 1 0.3 0.7 ];
- N = [ 2 1 ];
- [U R Q X] = qnsolve( "closed", N, QQ, V );
- assert( Q, [ 0.428 0.726 0.845; 0.108 0.158 0.734 ], 1e-3 );
- assert( X(1,1), 2.113, 1e-3 ); # CHECK
- assert( X(2,1), 0.524, 1e-3 ); # CHECK
- assert( all( all(U(:,[1,2])<=1) ) );
+ [r err] = dtmcchkP( [1 0 0; 0 0.5 0; 0 0 0] );
+ assert( r, 0 );
+ assert( index(err, "stochastic") > 0 );
 ***** test
- # example p. 26 Schwetman
- QQ = { qnmknode( "m/m/m-fcfs", [.25; .25] ),
-        qnmknode( "-/g/1-ps", [0; .1] ) };
- V = [1 0; 1 1];
- lambda = [1 0];
- N = [0 3];
- [U R Q X] = qnsolve( "mixed", lambda, N, QQ, V );
- assert( U(1,1), .25, 1e-3 );
- assert( X(1,1), 1.0, 1e-3 );
- assert( [R(1,1) R(2,1) R(2,2)], [1.201 0.885 0.135], 1e-3 );
-***** demo
- QQ = { qnmknode( "m/m/m-fcfs", [0.2 0.1 0.1; 0.2 0.1 0.1] ), ...
-        qnmknode( "-/g/1-ps", [0.4; 0.6] ), ...
-        qnmknode( "-/g/inf", [1; 2] ) };
- V = [ 1 0.6 0.4; ...
-       1 0.3 0.7 ];
- N = [ 2 1 ];
- [U R Q X] = qnsolve( "closed", N, QQ, V );
-9 tests, 9 passed, 0 known failure, 0 skipped
-[inst/qnmix.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnmix.m
+ P = [0 1; 1 0];
+ assert( dtmcchkP(P), 2 );
 ***** test
- lambda = [1 0 0];
- N = [1 1 1];
- S = V = [1 1 1; 1 1 1; 1 1 1];
- fail( "qnmix( lambda, N, S, V)", "same time");
- N = [0 0 1];
- fail( "qnmix( lambda, N, S, V)", "open nor closed" );
- N = [0 1 2];
- m = [ 1 1 2 ];
- fail( "qnmix( lambda, N, S, V, m)", "single-server and delay" );
- S = V = [1 1 1; 1 1 1];
- fail( "qnmix( lambda, N, S, V)", "rows" );
+ P = dtmcbd( linspace(0.1,0.4,10), linspace(0.4,0.1,10) );
+ assert( dtmcchkP(P), rows(P) );
 ***** test
- # Example p. 148 Zahorjan et al.
- lambda = [1 1/2 0 0];
- N = [0 0 1 1];
- V = [1 1; 1 1; 1 1; 1 1];
- S = [1/4 1/6; 1/2 1; 1/2 1; 1 4/3];
- [U R Q X] = qnmix(lambda, N, S, V );
- assert( Q(3,1), 4/19, 1e-4 );
- assert( Q(3,2), 15/19, 1e-4 );
- assert( Q(4,1), 5/19, 1e-4 );
- assert( Q(4,2), 14/19, 1e-4 );
- assert( Q, R.*X, 1e-5 ); # Little's Law
+ N = 1000;
+ P = reshape( 1:N^2, N, N );
+ P(1:N+1:end) = 0;
+ P = P ./ repmat(sum(P,2),1,N);
+ assert( dtmcchkP(P), N );
+5 tests, 5 passed, 0 known failure, 0 skipped
+[inst/dtmcisir.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/dtmcisir.m
 ***** test
- # Example 8.6 p. 345 Bolch et al.
- lambda = [0.5 0.25 0 0];
- N = [0 0 1 1];
- V = [2 1; 2.5 1.5; 1 0.5; 1 0.4];
- S = [0.4 0.6; 0.8 1.6; 0.3 0.5; 0.5 0.8];
- [U R Q X] = qnmix( lambda, N, S, V );
- assert( U([1 2],:), [0.4 0.3; 0.5 0.6], 1e-3 );
- assert( R([3 4],:), [4.829 6.951; 7.727 11.636], 1e-3 );
- assert( Q([3 4],:), [0.582 0.418; 0.624 0.376], 1e-3 );
- assert( Q([1 2],:), [8.822 5.383; 11.028 10.766], 1e-3 ); #FIXME
- assert( R([1 2],:), [8.822 10.766; 17.645 28.710], 1e-3 );
- assert( X(3,1)/V(3,1), 0.120, 1e-3 );
- assert( X(4,1)/V(4,1), 0.081, 1e-3 );
- assert( Q, R.*X, 1e-5 ); # Little's Law
+ P = [0 .5 0; 0 0 0];
+ fail( "dtmcisir(P)" );
 ***** test
- ## example figure 10 p. 26 Schwetman, "Implementing the Mean Value
- ## Analysis for the Solution of Queueing Network Models", Technical
- ## Report CSD-TR-355, feb 15, 1982, Purdue University. 
- S = [.25 0; .25 .10];
- V = [1 0; 1 1];
- lambda = [1 0];
- N = [0 3];
- [U R Q X] = qnmix( lambda, N, S, V );
- assert( U(1,1), .25, 1e-3 );
- assert( X(1,1), 1.0, 1e-3 );
- assert( [R(1,1) R(2,1) R(2,2)], [1.201 0.885 0.135], 1e-3 );
- assert( Q, R.*X, 1e-5 ); # Little's Law
-4 tests, 4 passed, 0 known failure, 0 skipped
-[inst/engset.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/engset.m
+ P = [0 1 0; 0 .5 .5; 0 1 0];
+ [r s] = dtmcisir(P);
+ assert( r == 0 );
+ assert( max(s), 2 );
+ assert( min(s), 1 );
 ***** test
- fail("erlangb(1, -1)", "positive");
- fail("erlangb(-1, 1)", "positive");
- fail("erlangb(1, 0)", "positive");
- fail("erlangb(0, 1)", "positive");
- fail("erlangb('foo',1)", "positive");
- fail("erlangb(1,'bar')", "positive");
- fail("erlangb([1 1],[1 1 1])","common size");
-1 test, 1 passed, 0 known failure, 0 skipped
+ P = [.5 .5 0; .2 .3 .5; 0 .2 .8];
+ [r s] = dtmcisir(P);
+ assert( r == 1 );
+ assert( max(s), 1 );
+ assert( min(s), 1 );
+3 tests, 3 passed, 0 known failure, 0 skipped
 [inst/ctmcchkQ.m]
 >>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/ctmcchkQ.m
 ***** test
@@ -6423,34 +6392,103 @@
  assert( result, 0 );
  assert( index(err, "infinitesimal") > 0 );
 5 tests, 5 passed, 0 known failure, 0 skipped
-[inst/qncspb.m]
->>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncspb.m
+[inst/qncscmva.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qncscmva.m
 ***** test
- fail( "qncspb( 1, [] )", "vector" );
- fail( "qncspb( 1, [0 -1])", "nonnegative" );
- fail( "qncspb( 0, [1 2] )", "positive" );
- fail( "qncspb( -1, [1 2])", "nonnegative" );
- fail( "qncspb( 1, [1 2], [1,1], [2, 2])", "single server" );
- fail( "qncspb( 1, [1 2], [1,1], [1, 1], -1)", "nonnegative" );
-***** # shared test function
-***** function test_pb( D, expected, Z=0 )
- for i=1:rows(expected)
-   N = expected(i,1);
-   [X_lower X_upper] = qncspb(N,D,ones(size(D)),ones(size(D)),Z);
-   X_exp_lower = expected(i,2);
-   X_exp_upper = expected(i,3);
-   assert( [N X_lower X_upper], [N X_exp_lower X_exp_upper], 1e-4 )
+ N=5;
+ S = [1 0.3 0.8 0.9];
+ V = [1 1 1 1];
+ [U1 R1 Q1 X1] = qncscmva( N, S(1:3), repmat(S(4),1,N), V );
+ [U2 R2 Q2 X2] = qncsmva(N, S, V);
+ assert( X1, X2, 1e-5 );
+ assert( U1, U2, 1e-5 );
+ assert( R1, R2, 1e-5 );
+ assert( Q1, Q2, 1e-5 );
+***** test
+ N=5;
+ S = [1 1 1 1 1; ...
+      1 1 1 1 1; ...
+      1 1 1 1 1; ...
+      1 1/2 1/3 1/4 1/5];
+ V = [1 1 1 1];
+ [U1 R1 Q1 X1] = qncscmva( N, S(1:3,1), S(4,:), V );
+ [U2 R2 Q2 X2] = qncsmvald(N, S, V);
+ assert( U1, U2, 1e-5 );
+ assert( R1, R2, 1e-5 );
+ assert( Q1, Q2, 1e-5 );
+ assert( X1, X2, 1e-5 );
+***** test
+ N=5;
+ S = [1 1 1 1 1; ...
+      1 1 1 1 1; ...
+      1 1 1 1 1; ...
+      1 1/2 1/3 1/4 1/5];
+ V = [1 2 1 1];
+ Z = 3;
+ [U1 R1 Q1 X1] = qncscmva( N, S(1:3,1), S(4,:), V, Z );
+ [U2 R2 Q2 X2] = qncsmvald(N, S, V, Z);
+ assert( U1, U2, 1e-5 );
+ assert( R1, R2, 1e-5 );
+ assert( Q1, Q2, 1e-5 );
+ assert( X1, X2, 1e-5 );
+***** demo
+ maxN = 90; # Max population size
+ Rmva = Rconv = Rcmva = zeros(1,maxN); # Results
+ S = 4; Z = 10; m = 8;
+ old = warning("query","qn:numerical-instability");
+ warning("off","qn:numerical-instability");
+ for N=1:maxN
+   [U R] = qncsmva(N,S,1,m,Z);		# Use MVA
+   Rmva(N) = R(1);
+   [U R] = qncsconv(N,[S Z],[1 1],[m -1]);	# Use Convolution
+   Rconv(N) = R(1);
+   if ( N > m )
+     Scmva = S ./ min(1:N,m);
+   else
+     Scmva = S ./ (1:N);
+   endif
+   [U R] = qncscmva(N,[],Scmva,1,Z);		# Use CMVA
+   Rcmva(N) = R(1);
  endfor
+ warning(old.state,"qn:numerical-instability");
+ plot(1:maxN, Rmva, ";MVA;", ...
+      1:maxN, Rconv, ";Convolution;", ...
+      1:maxN, Rcmva, ";CNVA;", "linewidth",2);
+ xlabel("Population size (N)");
+ ylabel("Response Time");
+ ax=axis(); ax(3) = 0; ax(4) = 40; axis(ax);
+ legend("location","northwest"); legend("boxoff");
+3 tests, 3 passed, 0 known failure, 0 skipped
+[inst/qnosvisits.m]
+>>>>> /build/reproducible-path/octave-queueing-1.2.7/inst/qnosvisits.m
 ***** test
- # table IV
- D = [ 0.1 0.1 0.09 0.08 ];
- #            N  X_lower  X_upper
- expected = [ 2  4.3174   4.3174; ... 
-              5  6.6600   6.7297; ...
-              10 8.0219   8.2700; ...
-              20 8.8672   9.3387; ...
-              80 9.6736   10.000 ];
- test_pb(D, expected);
+ fail( "qnosvisits([0 .5; .5 0],[0 -1])", "contains negative" );
+ fail( "qnosvisits([1 1 1; 1 1 1], [1 1])", "square" );
+***** test
+
+ ## Open, single class network
+
+ P = [0 0.2 0.5; 1 0 0; 1 0 0];
+ lambda = [ 0.1 0.3 0.2 ];
+ V = qnosvisits(P,lambda);
+ assert( V*P+lambda/sum(lambda),V,1e-5 );
+***** demo
+ p = 0.3;
+ lambda = 1.2
+ P = [0 0.3 0.5; ...
+      1 0   0  ; ...
+      1 0   0  ];
+ V = qnosvisits(P,[1.2 0 0])
+***** demo
+ P = [ 0 0.4 0.6 0; ...
+       0.2 0 0.2 0.6; ...
+       0 0 0 1; ...
+       0 0 0 0 ];
+ lambda = [0.1 0 0 0.3];
+ V = qnosvisits(P,lambda);
+ S = [2 1 2 1.8];
+ m = [3 1 1 2];
+ [U R Q X] = qnos( sum(lambda), S, V, m )
 2 tests, 2 passed, 0 known failure, 0 skipped
 Checking C++ files ...
 Done running the unit tests.
@@ -6481,12 +6519,14 @@
 dpkg-buildpackage: info: binary-only upload (no source included)
 dpkg-genchanges: info: not including original source code in upload
 I: copying local configuration
+I: user script /srv/workspace/pbuilder/11609/tmp/hooks/B01_cleanup starting
+I: user script /srv/workspace/pbuilder/11609/tmp/hooks/B01_cleanup finished
 I: unmounting dev/ptmx filesystem
 I: unmounting dev/pts filesystem
 I: unmounting dev/shm filesystem
 I: unmounting proc filesystem
 I: unmounting sys filesystem
 I: cleaning the build env 
-I: removing directory /srv/workspace/pbuilder/7751 and its subdirectories
-I: Current time: Sat Apr 27 19:57:52 -12 2024
-I: pbuilder-time-stamp: 1714291072
+I: removing directory /srv/workspace/pbuilder/11609 and its subdirectories
+I: Current time: Sun Apr 28 22:07:37 +14 2024
+I: pbuilder-time-stamp: 1714291657