with initial condition w(0) = v at
one or more time points.
Usage
## S4 method for signature 'matrix,vector'
expv(x, v, t = 1, u = NULL,
Markov = FALSE, transpose = Markov, ...)
## S4 method for signature 'CsparseMatrix,vector'
expv(x, v, t = 1,
u = NULL, Markov = FALSE, transpose = Markov, ...)
## S4 method for signature 'dgCMatrix,vector'
expv(x, v, t = 1, u = NULL,
Markov = FALSE, transpose = Markov, ...)
## S4 method for signature 'TsparseMatrix,vector'
expv(x, v, t = 1,
u = NULL, Markov = FALSE, transpose = Markov, ...)
## S4 method for signature 'matrix.csc,vector'
expv(x, v, t = 1, u = NULL,
Markov = FALSE, transpose = Markov, ...)
## S4 method for signature 'matrix.csr,vector'
expv(x, v, t = 1, u = NULL,
Markov = FALSE, transpose = Markov, ...)
## S4 method for signature 'matrix.coo,vector'
expv(x, v, t = 1, u = NULL,
Markov = FALSE, transpose = Markov, ...)
Arguments
x
a matrix.
v
a numeric vector. The initial value.
t
a numeric vector of time points at which
the solution is computed.
u
a numeric vector. Default NULL.
Markov
logical. If TRUE the matrix
is taken to be an rate matrix and steps are taken to
ensure that the computed result is a probability vector.
Default FALSE.
transpose
logical. If TRUE transpose
the matrix before the solution is computed. Default
equals Markov.
...
other arguments passed to
Rexpv.
Details
Analytically the solution is given as
w(t) = exp(tx)v + tphi(tx)u
with phi(z) =
(exp(z)-1)/z. For large matrices x the computation
of the full matrices exp(tx) and
phi(tx) is slow. An alternative is to
compute the solution w directly for a given initial
value v using Krylov subspace methods. This is, in
particular, efficient for large sparse matrices.
Note that if Q is a rate matrix for a homogeneous
continuous time Markov process (non-negative
off-diagonals and row-sums 0) and v is a
probability vector (the initial distribution at time
point 0), then the distribution at time point t
solves
w'(t) = Q^T w(t).
In this case we want to
take x to be Q^T to get the desired solution.
The solution is computed using the Fortran package
Expokit. Methods are available for matrix classes
implemented in the Matrix package as well as the SparseM
package. The implementation avoids the computation of the
full matrix exponential of tx and the approach is
advantageous when we want to compute w(t) for one
or a few initial values v. The full matrix
exponential should not be computed this way by
looping over n different initial values.
Though there is a method implemented for ordinary (dense)
matrices, such a matrix is simply coerced into a
CsparseMatrix before the solution is computed. It
is recommended that large sparse matrices are stored and
handled as such, e.g. using the classes and methods from
the Matrix package. Dense intermediates should be
avoided.
The x matrix is allowed to be a dense complex
matrix in which case v and u are also
allowed to be complex.
Value
An n by k matrix with k the length of
t and n the dimension of the matrix
x. The i'th column contains the solution of
the ODE at time point i.