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.63102 Finished Loading Parameters Finished ipoint() function Initial count = 4 upwindorderedv4.c:836:8: runtime error: pointer index expression with base 0x7fea1fd34800 overflowed to 0xf5f675e015c93df0 #0 0x7fea20103241 in updatetree /data/gannet/ripley/R/packages/tests-clang-SAN/QPot/src/upwindorderedv4.c:836:8 #1 0x7fea200fc704 in ordered_upwind /data/gannet/ripley/R/packages/tests-clang-SAN/QPot/src/upwindorderedv4.c:515:37 #2 0x7fea20105665 in quasipotential /data/gannet/ripley/R/packages/tests-clang-SAN/QPot/src/upwindorderedv4.c:1036:5 #3 0x768214 in do_dotCode /data/gannet/ripley/R/svn/R-devel/src/main/dotcode.c #4 0x86f932 in bcEval /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:6775:14 #5 0x853bb8 in Rf_eval /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:620:8 #6 0x8b0c5b in R_execClosure /data/gannet/ripley/R/svn/R-devel/src/main/eval.c #7 0x8adf35 in Rf_applyClosure /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:1706:16 #8 0x854ab8 in Rf_eval /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:743:12 #9 0x8bd1ab in do_set /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:2808:8 #10 0x85437b in Rf_eval /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:695:12 #11 0x988b56 in Rf_ReplIteration /data/gannet/ripley/R/svn/R-devel/src/main/main.c:260:2 #12 0x98d070 in R_ReplConsole /data/gannet/ripley/R/svn/R-devel/src/main/main.c:310:11 #13 0x98ce55 in run_Rmainloop /data/gannet/ripley/R/svn/R-devel/src/main/main.c:1108:5 #14 0x4da2ca in main /data/gannet/ripley/R/svn/R-devel/src/main/Rmain.c:29:5 #15 0x7fea2cc21f42 in __libc_start_main (/lib64/libc.so.6+0x23f42) #16 0x4302dd in _start (/data/gannet/ripley/R/R-clang-SAN/bin/exec/R+0x4302dd) SUMMARY: UndefinedBehaviorSanitizer: undefined-behavior upwindorderedv4.c:836:8 in upwindorderedv4.c:837:8: runtime error: pointer index expression with base 0x7fea1fd34800 overflowed to 0xf5f675e015c93df0 #0 0x7fea2010334c in updatetree /data/gannet/ripley/R/packages/tests-clang-SAN/QPot/src/upwindorderedv4.c:837:8 #1 0x7fea200fc704 in ordered_upwind /data/gannet/ripley/R/packages/tests-clang-SAN/QPot/src/upwindorderedv4.c:515:37 #2 0x7fea20105665 in quasipotential /data/gannet/ripley/R/packages/tests-clang-SAN/QPot/src/upwindorderedv4.c:1036:5 #3 0x768214 in do_dotCode /data/gannet/ripley/R/svn/R-devel/src/main/dotcode.c #4 0x86f932 in bcEval /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:6775:14 #5 0x853bb8 in Rf_eval /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:620:8 #6 0x8b0c5b in R_execClosure /data/gannet/ripley/R/svn/R-devel/src/main/eval.c #7 0x8adf35 in Rf_applyClosure /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:1706:16 #8 0x854ab8 in Rf_eval /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:743:12 #9 0x8bd1ab in do_set /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:2808:8 #10 0x85437b in Rf_eval /data/gannet/ripley/R/svn/R-devel/src/main/eval.c:695:12 #11 0x988b56 in Rf_ReplIteration /data/gannet/ripley/R/svn/R-devel/src/main/main.c:260:2 #12 0x98d070 in R_ReplConsole /data/gannet/ripley/R/svn/R-devel/src/main/main.c:310:11 #13 0x98ce55 in run_Rmainloop /data/gannet/ripley/R/svn/R-devel/src/main/main.c:1108:5 #14 0x4da2ca in main /data/gannet/ripley/R/svn/R-devel/src/main/Rmain.c:29:5 #15 0x7fea2cc21f42 in __libc_start_main (/lib64/libc.so.6+0x23f42) #16 0x4302dd in _start (/data/gannet/ripley/R/R-clang-SAN/bin/exec/R+0x4302dd) SUMMARY: UndefinedBehaviorSanitizer: undefined-behavior upwindorderedv4.c:837:8 in 4000 (2 31) is accepted, g=0.2514 Final count = 350 Finished ordered_upwind() function cputime = 0.97451 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.641925 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 = 0.971099 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.285882 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.431334 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.288255 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.440559 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.10291 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 = 2.58229 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.284557 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.427322 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.284162 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.435925 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.290292 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.4325 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.28212 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.432567 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.283526 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.428971 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.282834 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.445312 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.642209 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 = 0.952627 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.641683 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 = 0.957562 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) > > > > ### *