numint: Numerical integration of powers and fractions of powers via...

In episode: Estimation with Penalisation in Systems of Ordinary Differential Equations

Description

Evaluates numerical integrals of powers or fractions of powers of a d-dimensional function x.

Usage

numint(time, x, type, A, B)

Arguments

time	Vector (n) holding time points in between which the integrals are evaluated. Must be one series only (i.e., no decrements).
x	Matrix (mx(d+1)) holding discretised function. First column is time, must have no decrements. The remaining columns are function values, which will be interpolated.
type	String telling type, must be "power" or "fracpower".
A	Matrix (pxd) holding powers in rows.
B	Matrix (pxd) holding powers in rows.

Value

A matrix of dimension (m-1)xp holding the row-concatinated blocks of integrals: If type is "power" it evaluates the numerical integrals

(\int^{t_{i+1}}_{t_i}{x(s)^A})_i

where t_i are entries in first column of x and x(s) is constructed via a linear interpolation of x. If type is "fracpower" it evaluates the numerical integrals

(\int^{t_{i+1}}_{t_i}{x(s)^A / (1 + x(s)^B)})_i

where the fraction is evaluated coordinate wise.

The numerical integration uses the simpson rule using the intermediate time points in time in between any two consecutive time points in the first column of x. To get more accurate integrals include more intermediate time points in time.

See Also

imd, aim

Examples

# Trajectory of power law kinetics system
A <- matrix(c(1, 0, 1,
              1, 1, 0), byrow = TRUE, nrow = 2)
p <- plk(A)
x0 <- c(10, 4, 1)
theta <- matrix(c(0, -0.25,
                  0.75, 0,
                  0, -0.1), byrow = TRUE, nrow = 3)
Time <- seq(0, 1, by = .025)
traj <- numsolve(p, Time, x0, theta)

# Example: Integrate traj(s)^A between the time points in 'ti'
ti <- seq(0, 1, by = .1)
ss <- numint(time = ti, x = traj, type = "power", A = A, B = A)

# Example: Integrate traj(s)^A / (1 + traj(s)^B) between the time points in 'ti'
B <- matrix(c(0, 2, 1,
              2, 1, 0), byrow = TRUE, nrow = 2)
ss <- numint(time = ti, x = traj, type = "fracpower", A = A, B = B)

episode documentation built on May 1, 2019, 11:17 p.m.