Plot and Print Methods for Output of bvp solvers
Description
Plot the output of boundary value solver routines.
Usage
1 2 3 4 5 
Arguments
x 
the output of It is allowed to pass several objects of class

which 
the name(s) or the index to the variables that should be plotted. Default = all variables, except the first column. 
ask 
logical; if 
obs 
a By default the first column of an observed data set should contain
the If the first column of If 
obspar 
additional graphics arguments passed to 
... 
additional arguments. The graphical arguments are passed to The dots may also contain other objects of class
The arguments after ... must be matched exactly. 
Details
print.bvpSolve
prints the matrix without the attributes.
plot.bvpSolve
plots multiple figures on a page.
The number of panels per page is automatically determined up to 3 x 3
(par(mfrow = c(3, 3))
). This default can be overwritten by
specifying userdefined settings for mfrow
or mfcol
.
Set mfrow
equal to NULL
to avoid the plotting function to
change userdefined mfrow
or mfcol
settings.
Other graphical parameters can be passed as well. Parameters
are vectorized, either according to the number of plots
(xlab, ylab
, main, sub
, xlim, ylim
, log
,
asp, ann, axes, frame.plot
,panel.first,panel.last
,
cex.lab,cex.axis,cex.main
) or
according to the number of lines within one plot (other parameters
e.g. col
, lty
, lwd
etc.) so it is possible to
assign specific axis labels to individual plots, resp. different plotting
style. Plotting parameter ylim
, or xlim
can also be a list
to assign different axis limits to individual plots.
Similarly, the graphical parameters for observed data, as passed by
obspar
can be vectorized, according to the number of observed
data sets.
See Also
diagnostics.bvpSolve
, for a description of diagnostic messages.
Examples
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110  ## =============================================================================
## The example MUSN from Ascher et al., 1995.
## Numerical Solution of Boundary Value Problems for Ordinary Differential
## Equations, SIAM, Philadelphia, PA, 1995.
##
## The problem is
## u' = 0.5*u*(w  u)/v
## v' = 0.5*(w  u)
## w' = (0.9  1000*(w  y)  0.5*w*(w  u))/z
## z' = 0.5*(w  u)
## y' = 100*(y  w)
## on the interval [0 1] and with boundary conditions:
## u(0) = v(0) = w(0) = 1, z(0) = 10, w(1) = y(1)
## =============================================================================
musn < function(t, Y, pars) {
with (as.list(Y),
{
du < 0.5 * u * (wu)/v
dv < 0.5 * (wu)
dw < (0.9  1000 * (wy)  0.5 * w * (wu))/z
dz < 0.5 * (wu)
dy < 100 * (yw)
return(list(c(du, dv, dw, dz, dy)))
})
}
#
# Boundaries
#
bound < function(i,y,pars) {
with (as.list(y), {
if (i ==1) return (u1)
if (i ==2) return (v1)
if (i ==3) return (w1)
if (i ==4) return (z+10)
if (i ==5) return (wy)
})
}
xguess < seq(0, 1, len = 5)
yguess < matrix(ncol = 5, (rep(c(1, 1, 1, 10, 0.91), times = 5)) )
rownames(yguess) < c("u", "v", "w", "z", "y")
sol < bvpcol (bound = bound, x = seq(0, 1, by = 0.05),
leftbc = 4, func = musn, xguess = xguess, yguess = yguess)
mf < par("mfrow")
plot(sol)
par(mfrow = mf)
## =============================================================================
## Example 2. Example Problem 31 from Jeff Cash's website
## =============================================================================
Prob31 < function(t, Y, pars) {
with (as.list(Y), {
dy < sin(Tet)
dTet < M
dM < Q/xi
T < 1/cos (Tet) +xi*Q*tan(Tet)
dQ < 1/xi*((y1)*cos(Tet)M*T)
list(c( dy, dTet, dM, dQ))
})
}
ini < c(y = 0, Tet = NA, M = 0, Q = NA)
end < c(y = 0, Tet = NA, M = 0, Q = NA)
# run 1
xi <0.1
twp < bvptwp(yini = ini, yend = end, x = seq(0, 1, by = 0.01),
func = Prob31, atol = 1e10)
# run 2
xi < 0.05
twp2 < bvptwp(yini = ini, yend = end, x = seq(0, 1, by = 0.01),
func = Prob31, atol = 1e10)
# run 3
xi < 0.01
twp3 < bvptwp(yini = ini, yend = end, x = seq(0, 1, by = 0.01),
func = Prob31, atol = 1e10)
# print all outputs at once
plot(twp, twp2, twp3, xlab = "x", ylab = names(ini))
# change type, colors, ...
plot(twp, twp2, twp3, type = c("l", "b", "p"),
main = paste ("State Variable", names(ini)),
col = c("red", "blue", "orange"), cex = 2)
## =============================================================================
## Assume we have two 'data sets':
## =============================================================================
# data set in 'wide' format
obs1 < cbind(time = c(0, 0.5, 1), Tet = c(0.4, 0.0, 0.4))
# data set in 'long' format
obs2 < data.frame(name = "Tet", time = c(0, 0.5, 1), value = c(0.35, 0.0, 0.35))
plot(twp, twp2, obs = obs1, obspar = list(pch = 16, cex = 1.5))
plot(twp, twp2, obs = list(obs1, obs2),
obspar = list(pch = 16, cex = 1.5))
plot(twp, twp2, obs = list(obs1, obs2), which = c("Tet", "Q"),
obspar = list(pch = 16:17, cex = 1.5, col = c("red", "black"))
)
