ParallelModel {SoilR}

models for unconnected pools

Description

This function creates a (linear) numerical model for n independent (parallel) pools that can be queried afterwards. It is used by the convenient wrapper functions TwopParallelModel and ThreepParallelModel but can also be used independently.

Usage

ParallelModel(
  times,
  coeffs_tm,
  startvalues,
  inputrates,
  solverfunc = deSolve.lsoda.wrapper,
  pass = FALSE
)

Arguments

times

A vector containing the points in time where the solution is sought.

coeffs_tm

A TimeMap object consisting of a vector valued function containing the decay rates for the n pools as function of time and the time range where this function is valid. The length of the vector is equal to the number of pools.

startvalues

A vector containing the initial amount of carbon for the n pools. <<The length of this vector is equal to the number of pools and thus equal to the length of k. This is checked by the function.

inputrates

An object consisting of a vector valued function describing the inputs to the pools as functions of time TimeMap.new

solverfunc

The function used to actually solve the ODE system. This can be deSolve.lsoda.wrapper or any other user provided function with the same interface.

pass

if TRUE forces the constructor to create the model even if it is invalid

Examples

t_start=0 
t_end=10 
tn=50
timestep=(t_end-t_start)/tn 
t=seq(t_start,t_end,timestep) 
k=TimeMap(
function(times){c(-0.5,-0.2,-0.3)},
t_start,
t_end
)
c0=c(1, 2, 3)
#constant inputrates
inputrates=BoundInFluxes(
function(t){matrix(nrow=3,ncol=1,c(1,1,1))},
t_start,
t_end
) 
mod=ParallelModel(t,k,c0,inputrates)
Y=getC(mod)
lt1=1 ;lt2=2 ;lt3=3 
col1=1; col2=2; col3=3
plot(t,Y[,1],type="l",lty=lt1,col=col1,
ylab="C stocks",xlab="Time") 
lines(t,Y[,2],type="l",lty=lt2,col=col2) 
lines(t,Y[,3],type="l",lty=lt3,col=col3) 
legend(
"topleft",
c("C in pool 1",
"C in 2",
"C in pool 3"
),
lty=c(lt1,lt2,lt3),
col=c(col1,col2,col3)
)
Y=getAccumulatedRelease(mod)
plot(t,Y[,1],type="l",lty=lt1,col=col1,ylab="C release",xlab="Time") 
lines(t,Y[,2],lt2,type="l",lty=lt2,col=col2) 
lines(t,Y[,3],type="l",lty=lt3,col=col3) 
legend("topright",c("R1","R2","R3"),lty=c(lt1,lt2,lt3),col=c(col1,col2,col3))