R Under development (unstable) (2019-09-07 r77160) -- "Unsuffered Consequences" Copyright (C) 2019 The R Foundation for Statistical Computing Platform: x86_64-pc-linux-gnu (64-bit) R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. Natural language support but running in an English locale R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > pkgname <- "QPot" > source(file.path(R.home("share"), "R", "examples-header.R")) > options(warn = 1) > library('QPot') > > base::assign(".oldSearch", base::search(), pos = 'CheckExEnv') > base::assign(".old_wd", base::getwd(), pos = 'CheckExEnv') > cleanEx() > nameEx("Model2String") > ### * Model2String > > flush(stderr()); flush(stdout()) > > ### Name: Model2String > ### Title: Inserts parameter values into equations > ### Aliases: Model2String > > ### ** Examples > > # First example with the right hand side of an equation > test.eqn.x = "(alpha*x)*(1-(x/beta)) - ((delta*(x^2)*y)/(kappa + (x^2)))" > model.parms <- c(alpha=1.54, beta=10.14, delta=1, kappa=1) > equations.as.strings.x <- Model2String(test.eqn.x, parms = model.parms) Note: This function is supplied to help convert equations to strings. Long equations, equations spanning multiple lines, equations with strange notation, etc, may not work. Always check the output. [1] "(1.54*x)*(1-(x/10.14)) - ((1*(x^2)*y)/(1 + (x^2)))" > > # Second example with individual strings with left and right hand sides > # Note the use deSolve.form = TRUE to remove the dx and dy > test.eqn.x = "dx = (alpha*x)*(1-(x/beta)) - ((delta*(x^2)*y)/(kappa + (x^2)))" > test.eqn.y = "dy = ((gamma*(x^2)*y)/(kappa + (x^2))) - mu*(y^2)" > model.parms <- c(alpha=1.54, beta=10.14, delta=1, kappa=1, gamma=0.476, mu=0.112509) > equations.as.strings.x <- Model2String(test.eqn.x, parms = model.parms, + deSolve.form = TRUE, x.lhs.term = 'dx', y.lhs.term = 'dy') Note: This function is supplied to help convert equations to strings. Long equations, equations spanning multiple lines, equations with strange notation, etc, may not work. Always check the output. [1] "(1.54*x)*(1-(x/10.14)) - ((1*(x^2)*y)/(1 + (x^2)))" > equations.as.strings.y <- Model2String(test.eqn.y, parms = model.parms, + deSolve.form = TRUE, x.lhs.term = 'dx', y.lhs.term = 'dy') Note: This function is supplied to help convert equations to strings. Long equations, equations spanning multiple lines, equations with strange notation, etc, may not work. Always check the output. [1] "((0.476*(x^2)*y)/(1 + (x^2))) - 0.112509*(y^2)" > > # Third example with deSolve-style function call: > model.parms <- c(alpha=1.54, beta=10.14, delta=1, kappa=1, gamma=0.476, mu=0.112509) > ModelEquations <- function(t, state, parms) { + with(as.list(c(state, parms)), { + dx <- (alpha*x)*(1-(x/beta)) - ((delta*(x^2)*y)/(kappa + (x^2))) + dy <- ((gamma*(x^2)*y)/(kappa + (x^2))) - mu*(y^2) + list(c(dx,dy)) + }) + } > > Model2String(ModelEquations, parms = model.parms, deSolve.form = TRUE, + x.lhs.term = 'dx', y.lhs.term = 'dy') Note: This function is supplied to help convert equations to strings. Long equations, equations spanning multiple lines, equations with strange notation, etc, may not work. Always check the output. [1] " (1.54 * x) * (1 - (x/10.14)) - ((1 * (x^2) * y)/(1 + (x^2)))" [2] " ((0.476 * (x^2) * y)/(1 + (x^2))) - 0.112509 * (y^2)" [1] " (1.54 * x) * (1 - (x/10.14)) - ((1 * (x^2) * y)/(1 + (x^2)))" [2] " ((0.476 * (x^2) * y)/(1 + (x^2))) - 0.112509 * (y^2)" > > > > cleanEx() > nameEx("QPContour") > ### * QPContour > > flush(stderr()); flush(stdout()) > > ### Name: QPContour > ### Title: Contour plot of quasi-potential surfaces > ### Aliases: QPContour > > ### ** Examples > > # First, System of equations > equationx <- "1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x)" > equationy <- "((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y" > > # Second, shared parameters for each quasi-potential run > xbounds <- c(-0.5, 10.0) > ybounds <- c(-0.5, 10.0) > xstepnumber <- 150 > ystepnumber <- 150 > > # Third, first local quasi-potential run > xinit1 <- 1.40491 > yinit1 <- 2.80808 > storage.eq1 <- QPotential(x.rhs = equationx, x.start = xinit1, + x.bound = xbounds, x.num.steps = xstepnumber, y.rhs = equationy, + y.start = yinit1, y.bound = ybounds, y.num.steps = ystepnumber) The upwind ordered method will be chatty Set verboseC = FALSE to turn off completely Creating file name. File name created. Completed Memory Allocation equationx = 1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x) equationy = ((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y hx = 0.0704698 hy = 0.0704698 Finished initializing a bunch of matrices in param() function cputime = 0.673692 Finished Loading Parameters Finished ipoint() function Initial count = 4 upwindorderedv4.c:836:9: runtime error: pointer index expression with base 0x7f12e8fc5800 overflowed to 0xf5f67508def24df0 #0 0x7f12e935a00f in updatetree /data/gannet/ripley/R/packages/tests-gcc-SAN/QPot/src/upwindorderedv4.c:836 #1 0x7f12e935d130 in ordered_upwind /data/gannet/ripley/R/packages/tests-gcc-SAN/QPot/src/upwindorderedv4.c:515 #2 0x7f12e935fb83 in quasipotential /data/gannet/ripley/R/packages/tests-gcc-SAN/QPot/src/upwindorderedv4.c:1036 #3 0x5845c3 in do_dotCode /data/gannet/ripley/R/svn/R-devel/src/main/dotcode.c:1854 #4 0x612077 in bcEval /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:6775 #5 0x655017 in Rf_eval /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:620 #6 0x65a5c5 in R_execClosure /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:1780 #7 0x65ca08 in Rf_applyClosure /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:1706 #8 0x655448 in Rf_eval /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:743 #9 0x6616c1 in do_set /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:2808 #10 0x65597c in Rf_eval /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:695 #11 0x6ce625 in Rf_ReplIteration /data/gannet/ripley/R/svn/R-devel/src/main/main.c:260 #12 0x6ce625 in Rf_ReplIteration /data/gannet/ripley/R/svn/R-devel/src/main/main.c:200 #13 0x6ced18 in R_ReplConsole /data/gannet/ripley/R/svn/R-devel/src/main/main.c:310 #14 0x6cee64 in run_Rmainloop /data/gannet/ripley/R/svn/R-devel/src/main/main.c:1108 #15 0x4192e8 in main /data/gannet/ripley/R/svn/R-devel/src/main/Rmain.c:29 #16 0x7f12f5752f42 in __libc_start_main (/lib64/libc.so.6+0x23f42) #17 0x41ba1d in _start (/data/gannet/ripley/R/gcc-SAN/bin/exec/R+0x41ba1d) upwindorderedv4.c:837:9: runtime error: pointer index expression with base 0x7f12e8fc5800 overflowed to 0xf5f67508def24df0 #0 0x7f12e9359f93 in updatetree /data/gannet/ripley/R/packages/tests-gcc-SAN/QPot/src/upwindorderedv4.c:837 #1 0x7f12e935d130 in ordered_upwind /data/gannet/ripley/R/packages/tests-gcc-SAN/QPot/src/upwindorderedv4.c:515 #2 0x7f12e935fb83 in quasipotential /data/gannet/ripley/R/packages/tests-gcc-SAN/QPot/src/upwindorderedv4.c:1036 #3 0x5845c3 in do_dotCode /data/gannet/ripley/R/svn/R-devel/src/main/dotcode.c:1854 #4 0x612077 in bcEval /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:6775 #5 0x655017 in Rf_eval /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:620 #6 0x65a5c5 in R_execClosure /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:1780 #7 0x65ca08 in Rf_applyClosure /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:1706 #8 0x655448 in Rf_eval /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:743 #9 0x6616c1 in do_set /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:2808 #10 0x65597c in Rf_eval /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:695 #11 0x6ce625 in Rf_ReplIteration /data/gannet/ripley/R/svn/R-devel/src/main/main.c:260 #12 0x6ce625 in Rf_ReplIteration /data/gannet/ripley/R/svn/R-devel/src/main/main.c:200 #13 0x6ced18 in R_ReplConsole /data/gannet/ripley/R/svn/R-devel/src/main/main.c:310 #14 0x6cee64 in run_Rmainloop /data/gannet/ripley/R/svn/R-devel/src/main/main.c:1108 #15 0x4192e8 in main /data/gannet/ripley/R/svn/R-devel/src/main/Rmain.c:29 #16 0x7f12f5752f42 in __libc_start_main (/lib64/libc.so.6+0x23f42) #17 0x41ba1d in _start (/data/gannet/ripley/R/gcc-SAN/bin/exec/R+0x41ba1d) 4000 (2 31) is accepted, g=0.2514 Final count = 350 Finished ordered_upwind() function cputime = 1.26021 Saves only to R In datasave case 2 Successful. Exiting C code > > # Fourth, second local quasi-potential run > xinit2 <- 4.9040 > yinit2 <- 4.06187 > storage.eq2 <- QPotential(x.rhs = equationx, x.start = xinit2, + x.bound = xbounds, x.num.steps = xstepnumber, y.rhs = equationy, + y.start = yinit2, y.bound = ybounds, y.num.steps = ystepnumber) The upwind ordered method will be chatty Set verboseC = FALSE to turn off completely Creating file name. File name created. Completed Memory Allocation equationx = 1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x) equationy = ((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y hx = 0.0704698 hy = 0.0704698 Finished initializing a bunch of matrices in param() function cputime = 0.650772 Finished Loading Parameters Finished ipoint() function Initial count = 4 4003 (2 31) is accepted, g=0.2340 Final count = 350 Finished ordered_upwind() function cputime = 1.22894 Saves only to R In datasave case 2 Successful. Exiting C code > > # Fifth, determine global quasi-potential > unst.x <- c(0, 4.2008) > unst.y <- c(0, 4.0039) > ex1.global <- QPGlobal(local.surfaces = list(storage.eq1, storage.eq2), + unstable.eq.x = unst.x, unstable.eq.y = unst.y, x.bound = xbounds, + y.bound = ybounds) > > # Sixth, contour of the quasi-potential > QPContour(ex1.global, dens = c(100,100), x.bound = xbounds, + y.bound = ybounds, c.parm = 5) > > > > cleanEx() > nameEx("QPGlobal") > ### * QPGlobal > > flush(stderr()); flush(stdout()) > > ### Name: QPGlobal > ### Title: Finding the global quasi-potential > ### Aliases: QPGlobal > ### Keywords: Global quasi-potential > > ### ** Examples > > # First, System of equations > equationx <- "1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x)" > equationy <- "((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y" > > # Second, shared parameters for each quasi-potential run > xbounds <- c(-0.5, 10.0) > ybounds <- c(-0.5, 10.0) > xstepnumber <- 100 > ystepnumber <- 100 > > # Third, first local quasi-potential run > xinit1 <- 1.40491 > yinit1 <- 2.80808 > storage.eq1 <- QPotential(x.rhs = equationx, x.start = xinit1, + x.bound = xbounds, x.num.steps = xstepnumber, y.rhs = equationy, + y.start = yinit1, y.bound = ybounds, y.num.steps = ystepnumber) The upwind ordered method will be chatty Set verboseC = FALSE to turn off completely Creating file name. File name created. Completed Memory Allocation equationx = 1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x) equationy = ((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y hx = 0.106061 hy = 0.106061 Finished initializing a bunch of matrices in param() function cputime = 0.27894 Finished Loading Parameters Finished ipoint() function Initial count = 4 1722 (2 15) is accepted, g=0.2382 Final count = 246 Finished ordered_upwind() function cputime = 0.544095 Saves only to R In datasave case 2 Successful. Exiting C code > > # Fourth, second local quasi-potential run > xinit2 <- 4.9040 > yinit2 <- 4.06187 > storage.eq2 <- QPotential(x.rhs = equationx, x.start = xinit2, + x.bound = xbounds, x.num.steps = xstepnumber, y.rhs = equationy, + y.start = yinit2, y.bound = ybounds, y.num.steps = ystepnumber) The upwind ordered method will be chatty Set verboseC = FALSE to turn off completely Creating file name. File name created. Completed Memory Allocation equationx = 1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x) equationy = ((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y hx = 0.106061 hy = 0.106061 Finished initializing a bunch of matrices in param() function cputime = 0.285765 Finished Loading Parameters Finished ipoint() function Initial count = 4 1721 (2 15) is accepted, g=0.2236 Final count = 246 Finished ordered_upwind() function cputime = 0.550071 Saves only to R In datasave case 2 Successful. Exiting C code > > # Fifth, determine global quasi-potential > unst.x <- c(0, 4.2008) > unst.y <- c(0, 4.0039) > ex1.global <- QPGlobal(local.surfaces = list(storage.eq1, storage.eq2), + unstable.eq.x = unst.x, unstable.eq.y = unst.y, x.bound = xbounds, + y.bound = ybounds) > > > > cleanEx() > nameEx("QPotential") > ### * QPotential > > flush(stderr()); flush(stdout()) > > ### Name: QPotential > ### Title: Computes the quasi-potential for a system of stochastic > ### differential equations using the ordered upwind method. > ### Aliases: QPotential > > ### ** Examples > > # First, System of equations > equationx <- "1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x)" > equationy <- "((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y" > > # Second, shared parameters for each quasi-potential run > xbounds <- c(-0.5, 8.0) > ybounds <- c(-0.5, 8.0) > xstepnumber <- 200 > ystepnumber <- 200 > > # Third, a local quasi-potential run > xinit1 <- 1.40491 > yinit1 <- 2.80808 > storage.eq1 <- QPotential(x.rhs = equationx, x.start = xinit1, + x.bound = xbounds, x.num.steps = xstepnumber, y.rhs = equationy, + y.start = yinit1, y.bound = ybounds, y.num.steps = ystepnumber) The upwind ordered method will be chatty Set verboseC = FALSE to turn off completely Creating file name. File name created. Completed Memory Allocation equationx = 1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x) equationy = ((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y hx = 0.0427136 hy = 0.0427136 Finished initializing a bunch of matrices in param() function cputime = 1.14789 Finished Loading Parameters Finished ipoint() function Initial count = 4 11538 (2 49) is accepted, g=0.2872 Final count = 588 Finished ordered_upwind() function cputime = 3.15321 Saves only to R In datasave case 2 Successful. Exiting C code > # Visualize the quasi-potential > QPContour(storage.eq1, dens = c(xstepnumber, ystepnumber), + x.bound = xbounds, y.bound = ybounds, c.parm = 5) > > > > cleanEx() > nameEx("TSDensity") > ### * TSDensity > > flush(stderr()); flush(stdout()) > > ### Name: TSDensity > ### Title: Density plot from simulation of two-dimensional stochastic > ### differential equations > ### Aliases: TSDensity > ### Keywords: plot simulations stochastic > > ### ** Examples > > ## Not run: > ##D # First, the parameter values, as found in TSTraj > ##D model.state <- c(x = 3, y = 3) > ##D model.sigma <- 0.2 > ##D model.deltat <- 0.005 > ##D model.time <- 100 > ##D > ##D # Second, write out the deterministic skeleton of the equations to be simulated, > ##D # as found in TSTraj > ##D #Example 1 from article > ##D equationx <- "1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0 + x*x)" > ##D equationy <- "((0.476*x*x*y)/(1 + x*x)) - 0.112590*y*y" > ##D > ##D # Third, run it, as found in TSTraj > ##D ModelOut <- TSTraj(y0 = model.state, time = model.time, deltat = model.deltat, > ##D x.rhs = equationx, y.rhs = equationy, sigma = model.sigma) > ##D # Fourth, plot it: > ##D # in 1D > ##D TSDensity(ModelOut, dim = 1) > ##D # in 2D > ##D TSDensity(ModelOut, dim = 2, kde2d.n = 20, xlab = "") > ## End(Not run) > > > > cleanEx() > nameEx("TSPlot") > ### * TSPlot > > flush(stderr()); flush(stdout()) > > ### Name: TSPlot > ### Title: Plot simulation of two-dimensional stochastic differential > ### equations > ### Aliases: TSPlot > ### Keywords: plot simulations stochastic > > ### ** Examples > > # First, the parameter values, as found in TSTraj > model.state <- c(x = 3, y = 3) > model.sigma <- 0.2 > model.deltat <- 0.05 > model.time <- 100 > > # Second, write out the deterministic skeleton of the equations to be simulated, > # as found in TSTraj > #Example 1 from article > equationx <- "1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0 + x*x)" > equationy <- "((0.476*x*x*y)/(1 + x*x)) - 0.112590*y*y" > > # Third, run it, as found in TSTraj > ModelOut <- TSTraj(y0 = model.state, time = model.time, deltat = model.deltat, + x.rhs = equationx, y.rhs = equationy, sigma = model.sigma) > # Fourth, plot it: > # in 1D > TSPlot(ModelOut, deltat = model.deltat, dim = 1) > # in 2D > TSPlot(ModelOut, deltat = model.deltat, dim = 2) > > > > cleanEx() > nameEx("TSTraj") > ### * TSTraj > > flush(stderr()); flush(stdout()) > > ### Name: TSTraj > ### Title: Simulate two-dimensional stochastic differential equations > ### Aliases: TSTraj > ### Keywords: Stochastic simulation > > ### ** Examples > > # First, the parameter values > model.state <- c(x = 3, y = 3) > model.sigma <- 0.2 > model.deltat <- 0.1 > model.time <- 100 > > # Second, write out the deterministic skeleton of the equations to be simulated > equationx <- "1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0 + x*x)" > equationy <- "((0.476*x*x*y)/(1 + x*x)) - 0.112590*y*y" > > # Third, Run it > ModelOut <- TSTraj(y0 = model.state, time = model.time, deltat = model.deltat, + x.rhs = equationx, y.rhs = equationy, sigma = model.sigma) > > # Can also input x.rhs and y.rhs as strings that contain parameter names > # and include parms with names and values of parameters > model.state <- c(x = 1, y = 2) > model.parms <- c(alpha = 1.54, beta = 10.14, delta = 1, kappa = 1, gamma = 0.476, mu = 0.112509) > model.sigma <- 0.2 > model.time <- 100 > model.deltat <- 0.1 > > test.eqn.x = "(alpha*x)*(1-(x/beta)) - ((delta*(x^2)*y)/(kappa + (x^2)))" > test.eqn.y = "((gamma*(x^2)*y)/(kappa + (x^2))) - mu*(y^2)" > > ModelOut.parms <- TSTraj(y0 = model.state, time = model.time, deltat = model.deltat, + x.rhs = test.eqn.x, y.rhs = test.eqn.y, parms = model.parms, sigma = model.sigma) > > > > cleanEx() > nameEx("VecDecomAll") > ### * VecDecomAll > > flush(stderr()); flush(stdout()) > > ### Name: VecDecomAll > ### Title: Vector decomposition and remainder fields > ### Aliases: VecDecomAll > > ### ** Examples > > # First, the system of equations > equationx <- "1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x)" > equationy <- "((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y" > > # Second, shared parameters for each quasi-potential run > xbounds <- c(-0.5, 10.0) > ybounds <- c(-0.5, 10.0) > xstepnumber <- 100 > ystepnumber <- 100 > > # Third, first local quasi-potential run > xinit1 <- 1.40491 > yinit1 <- 2.80808 > storage.eq1 <- QPotential(x.rhs = equationx, x.start = xinit1, + x.bound = xbounds, x.num.steps = xstepnumber, y.rhs = equationy, + y.start = yinit1, y.bound = ybounds, y.num.steps = ystepnumber) The upwind ordered method will be chatty Set verboseC = FALSE to turn off completely Creating file name. File name created. Completed Memory Allocation equationx = 1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x) equationy = ((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y hx = 0.106061 hy = 0.106061 Finished initializing a bunch of matrices in param() function cputime = 0.282295 Finished Loading Parameters Finished ipoint() function Initial count = 4 1722 (2 15) is accepted, g=0.2382 Final count = 246 Finished ordered_upwind() function cputime = 0.542937 Saves only to R In datasave case 2 Successful. Exiting C code > > # Fourth, second local quasi-potential run > xinit2 <- 4.9040 > yinit2 <- 4.06187 > storage.eq2 <- QPotential(x.rhs = equationx, x.start = xinit2, + x.bound = xbounds, x.num.steps = xstepnumber, y.rhs = equationy, + y.start = yinit2, y.bound = ybounds, y.num.steps = ystepnumber) The upwind ordered method will be chatty Set verboseC = FALSE to turn off completely Creating file name. File name created. Completed Memory Allocation equationx = 1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x) equationy = ((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y hx = 0.106061 hy = 0.106061 Finished initializing a bunch of matrices in param() function cputime = 0.295936 Finished Loading Parameters Finished ipoint() function Initial count = 4 1721 (2 15) is accepted, g=0.2236 Final count = 246 Finished ordered_upwind() function cputime = 0.551742 Saves only to R In datasave case 2 Successful. Exiting C code > > # Fifth, determine global quasi-potential > unst.x <- c(0, 4.2008) > unst.y <- c(0, 4.0039) > ex1.global <- QPGlobal(local.surfaces = list(storage.eq1, storage.eq2), + unstable.eq.x = unst.x, unstable.eq.y = unst.y, x.bound = xbounds, + y.bound = ybounds) > > # Sixth, decompose the global quasi-potential into the > # deterministic skeleton, gradient, and remainder vector fields > VDAll <- VecDecomAll(surface = ex1.global, x.rhs = equationx, y.rhs = equationy, + x.bound = xbounds, y.bound = ybounds) > > > > cleanEx() > nameEx("VecDecomGrad") > ### * VecDecomGrad > > flush(stderr()); flush(stdout()) > > ### Name: VecDecomGrad > ### Title: Vector decomposition and remainder fields > ### Aliases: VecDecomGrad > ### Keywords: decomposition, field gradient vector > > ### ** Examples > > # First, the system of equations > equationx <- "1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x)" > equationy <- "((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y" > > # Second, shared parameters for each quasi-potential run > xbounds <- c(-0.5, 10.0) > ybounds <- c(-0.5, 10.0) > xstepnumber <- 100 > ystepnumber <- 100 > > # Third, first local quasi-potential run > xinit1 <- 1.40491 > yinit1 <- 2.80808 > storage.eq1 <- QPotential(x.rhs = equationx, x.start = xinit1, + x.bound = xbounds, x.num.steps = xstepnumber, y.rhs = equationy, + y.start = yinit1, y.bound = ybounds, y.num.steps = ystepnumber) The upwind ordered method will be chatty Set verboseC = FALSE to turn off completely Creating file name. File name created. Completed Memory Allocation equationx = 1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x) equationy = ((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y hx = 0.106061 hy = 0.106061 Finished initializing a bunch of matrices in param() function cputime = 0.295321 Finished Loading Parameters Finished ipoint() function Initial count = 4 1722 (2 15) is accepted, g=0.2382 Final count = 246 Finished ordered_upwind() function cputime = 0.546363 Saves only to R In datasave case 2 Successful. Exiting C code > > # Fourth, second local quasi-potential run > xinit2 <- 4.9040 > yinit2 <- 4.06187 > storage.eq2 <- QPotential(x.rhs = equationx, x.start = xinit2, + x.bound = xbounds, x.num.steps = xstepnumber, y.rhs = equationy, + y.start = yinit2, y.bound = ybounds, y.num.steps = ystepnumber) The upwind ordered method will be chatty Set verboseC = FALSE to turn off completely Creating file name. File name created. Completed Memory Allocation equationx = 1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x) equationy = ((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y hx = 0.106061 hy = 0.106061 Finished initializing a bunch of matrices in param() function cputime = 0.271439 Finished Loading Parameters Finished ipoint() function Initial count = 4 1721 (2 15) is accepted, g=0.2236 Final count = 246 Finished ordered_upwind() function cputime = 0.555925 Saves only to R In datasave case 2 Successful. Exiting C code > > # Fifth, determine global quasi-potential > unst.x <- c(0, 4.2008) > unst.y <- c(0, 4.0039) > ex1.global <- QPGlobal(local.surfaces = list(storage.eq1, storage.eq2), + unstable.eq.x = unst.x, unstable.eq.y = unst.y, x.bound = xbounds, + y.bound = ybounds) > > # Sixth, creat the gradient vector field > VDG <- VecDecomGrad(surface = ex1.global) > > > > cleanEx() > nameEx("VecDecomPlot") > ### * VecDecomPlot > > flush(stderr()); flush(stdout()) > > ### Name: VecDecomPlot > ### Title: Plotting function for vector decomposition and remainder fields > ### Aliases: VecDecomPlot > ### Keywords: detrministic field gradient plot plot, remainder skeleton > ### vector > > ### ** Examples > > # First, system of equations > equationx <- "1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x)" > equationy <- "((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y" > > # Second, shared parameters for each quasi-potential run > xbounds <- c(-0.5, 10.0) > ybounds <- c(-0.5, 10.0) > xstepnumber <- 100 > ystepnumber <- 100 > > # Third, first local quasi-potential run > xinit1 <- 1.40491 > yinit1 <- 2.80808 > storage.eq1 <- QPotential(x.rhs = equationx, x.start = xinit1, + x.bound = xbounds, x.num.steps = xstepnumber, y.rhs = equationy, + y.start = yinit1, y.bound = ybounds, y.num.steps = ystepnumber) The upwind ordered method will be chatty Set verboseC = FALSE to turn off completely Creating file name. File name created. Completed Memory Allocation equationx = 1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x) equationy = ((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y hx = 0.106061 hy = 0.106061 Finished initializing a bunch of matrices in param() function cputime = 0.304514 Finished Loading Parameters Finished ipoint() function Initial count = 4 1722 (2 15) is accepted, g=0.2382 Final count = 246 Finished ordered_upwind() function cputime = 0.542919 Saves only to R In datasave case 2 Successful. Exiting C code > > # Fourth, second local quasi-potential run > xinit2 <- 4.9040 > yinit2 <- 4.06187 > storage.eq2 <- QPotential(x.rhs = equationx, x.start = xinit2, + x.bound = xbounds, x.num.steps = xstepnumber, y.rhs = equationy, + y.start = yinit2, y.bound = ybounds, y.num.steps = ystepnumber) The upwind ordered method will be chatty Set verboseC = FALSE to turn off completely Creating file name. File name created. Completed Memory Allocation equationx = 1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x) equationy = ((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y hx = 0.106061 hy = 0.106061 Finished initializing a bunch of matrices in param() function cputime = 0.288241 Finished Loading Parameters Finished ipoint() function Initial count = 4 1721 (2 15) is accepted, g=0.2236 Final count = 246 Finished ordered_upwind() function cputime = 0.571706 Saves only to R In datasave case 2 Successful. Exiting C code > > # Fifth, determine global quasi-potential > unst.x <- c(0, 4.2008) > unst.y <- c(0, 4.0039) > ex1.global <- QPGlobal(local.surfaces = list(storage.eq1, storage.eq2), + unstable.eq.x = unst.x, unstable.eq.y = unst.y, x.bound = xbounds, + y.bound = ybounds) > > # Sixth, decompose the global quasi-potential into the > # deterministic skeleton, gradient, and remainder vector fields > VDAll <- VecDecomAll(surface = ex1.global, x.rhs = equationx, y.rhs = equationy, + x.bound = xbounds, y.bound = ybounds) > > # Seventh, plot all three vector fields > # The deterministic skeleton vector field > VecDecomPlot(x.field = VDAll[,,1], y.field = VDAll[,,2], dens = c(25,25), + x.bound = xbounds, y.bound = ybounds, tail.length = 0.25, head.length = 0.05) > # The gradient vector field > VecDecomPlot(x.field = VDAll[,,3], y.field = VDAll[,,4], dens = c(25,25), + x.bound = xbounds, y.bound = ybounds, tail.length = 0.15, head.length = 0.05) > # The remainder vector field > VecDecomPlot(x.field = VDAll[,,5], y.field = VDAll[,,6], dens = c(25,25), + x.bound = xbounds, y.bound = ybounds, tail.length = 0.15, head.length = 0.05) > > > > cleanEx() > nameEx("VecDecomRem") > ### * VecDecomRem > > flush(stderr()); flush(stdout()) > > ### Name: VecDecomRem > ### Title: Vector decomposition and remainder fields > ### Aliases: VecDecomRem > ### Keywords: decompoosition, field remainder vector > > ### ** Examples > > # First, the system of equations > equationx <- "1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x)" > equationy <- "((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y" > > # Second, shared parameters for each quasi-potential run > xbounds <- c(-0.5, 10.0) > ybounds <- c(-0.5, 10.0) > xstepnumber <- 150 > ystepnumber <- 150 > > # Third, first local quasi-potential run > xinit1 <- 1.40491 > yinit1 <- 2.80808 > storage.eq1 <- QPotential(x.rhs = equationx, x.start = xinit1, + x.bound = xbounds, x.num.steps = xstepnumber, y.rhs = equationy, + y.start = yinit1, y.bound = ybounds, y.num.steps = ystepnumber) The upwind ordered method will be chatty Set verboseC = FALSE to turn off completely Creating file name. File name created. Completed Memory Allocation equationx = 1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x) equationy = ((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y hx = 0.0704698 hy = 0.0704698 Finished initializing a bunch of matrices in param() function cputime = 0.656378 Finished Loading Parameters Finished ipoint() function Initial count = 4 4000 (2 31) is accepted, g=0.2514 Final count = 350 Finished ordered_upwind() function cputime = 1.18682 Saves only to R In datasave case 2 Successful. Exiting C code > > # Fourth, second local quasi-potential run > xinit2 <- 4.9040 > yinit2 <- 4.06187 > storage.eq2 <- QPotential(x.rhs = equationx, x.start = xinit2, + x.bound = xbounds, x.num.steps = xstepnumber, y.rhs = equationy, + y.start = yinit2, y.bound = ybounds, y.num.steps = ystepnumber) The upwind ordered method will be chatty Set verboseC = FALSE to turn off completely Creating file name. File name created. Completed Memory Allocation equationx = 1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x) equationy = ((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y hx = 0.0704698 hy = 0.0704698 Finished initializing a bunch of matrices in param() function cputime = 0.65172 Finished Loading Parameters Finished ipoint() function Initial count = 4 4003 (2 31) is accepted, g=0.2340 Final count = 350 Finished ordered_upwind() function cputime = 1.22387 Saves only to R In datasave case 2 Successful. Exiting C code > > # Fifth, determine global quasi-potential > unst.x <- c(0, 4.2008) > unst.y <- c(0, 4.0039) > ex1.global <- QPGlobal(local.surfaces = list(storage.eq1, storage.eq2), + unstable.eq.x = unst.x, unstable.eq.y = unst.y, x.bound = xbounds, + y.bound = ybounds) > > # Sixth, create remainder field > VDR <- VecDecomRem(surface = ex1.global, x.rhs = equationx, y.rhs = equationy, + x.bound = xbounds, y.bound = ybounds) > > > > cleanEx() > nameEx("VecDecomVec") > ### * VecDecomVec > > flush(stderr()); flush(stdout()) > > ### Name: VecDecomVec > ### Title: Vector decomposition and remainder fields > ### Aliases: VecDecomVec > > ### ** Examples > > # First, the system of equations > equationx <- "1.54*x*(1.0-(x/10.14)) - (y*x*x)/(1.0+x*x)" > equationy <- "((0.476*x*x*y)/(1+x*x)) - 0.112590*y*y" > > # Second, shared parameters for each quasi-potential run > xbounds <- c(-0.5, 20.0) > ybounds <- c(-0.5, 20.0) > xstepnumber <- 1000 > ystepnumber <- 1000 > > # Third, create the deterministic skeleton vector field > VDV <- VecDecomVec(x.num.steps = xstepnumber, y.num.steps = ystepnumber, x.rhs = equationx, + y.rhs = equationy, x.bound = xbounds, y.bound = ybounds) > > > > ### *