ctmc_kms.h File Reference

Koury-McAllister-Stewart aggregation-disaggregation for a nearly completely decomposable CTMC. More...

#include <cstddef>
#include <vector>
#include "line/api/mc/ctmc_courtois.h"
#include "line/api/mc/ctmc_bicgstab.h"
#include "line/api/mc/ctmc_gmres.h"
#include "line/api/mc/dtmc_solve.h"
#include "line/num/number.h"
#include "line/util/error.h"
#include "line/util/lu.h"
#include "line/util/matrix.h"

Include dependency graph for ctmc_kms.h:

Go to the source code of this file.

Classes
struct	line::mc::KmsResult< T >

Namespaces
namespace	line
namespace	line::mc

Functions
template<class T>
KmsResult< T >	line::mc::ctmc_kms (const Matrix< T > &Q, const std::vector< std::vector< std::size_t > > &MS, std::size_t numSteps)
	Koury-McAllister-Stewart aggregation-disaggregation for a nearly completely decomposable CTMC.

Detailed Description

Koury-McAllister-Stewart aggregation-disaggregation for a nearly completely decomposable CTMC.

Templated port of matlab/src/api/mc/ctmc_kms.m and jar/src/main/java/jline/api/mc/Ctmc_kms.java. Starting from the Courtois approximation, each sweep conditions the current iterate within every macro-state, aggregates those conditional vectors into a macro-state chain, solves it for the macro-state weights, and disaggregates by one block Gauss-Seidel sweep on the uniformized matrix:

pn (D - U) = zn L, D = blkdiag(I - P), L and U the strictly block lower and upper parts of P, zn the reweighted conditional.

Everything is done in the PERMUTED, macro-state-major index space, because that is the space P and the block offsets live in; ctmc_courtois returns its vector already mapped back to the original ordering, so the initial iterate is permuted on entry and the result is unpermuted on exit. Mixing the two spaces is the bug this port is written to avoid: with contiguous macro-states the two orderings coincide and the error is invisible.

As in the reference there is no convergence test: the caller asks for a fixed number of sweeps. Above 6000 states the block solve is attempted by GMRES first, exactly as MATLAB does, since that is where the direct factorization stops fitting in memory; a non-zero GMRES flag falls back to the direct solve, so the result never depends on whether the iteration converged.

GATED ON TRANSCENDENTAL ARITHMETIC: it is seeded by ctmc_courtois, whose epsMAX is an eigenvalue modulus, and its own large-system path is GMRES.

Definition in file ctmc_kms.h.