[R-sig-dyn-mod] Problem with simecol and ode - crash from R

Soetaert, Karline K.Soetaert at nioo.knaw.nl
Wed Oct 21 18:36:26 CEST 2009


Rainer,
 
I cannot recreate your crash; I get a lot of warnings and finally an error from lsoda, which terminates the simulation (but not R).  
The warnings are there because, between 43 and 44 simulation units the state variable becomes -Inf.  This is due to the 0-order consumption rate (h), so that "N" continues being consumed even if there is nothing to consume. This generates negative densities, which are also being consumed etc..., and N decreases exponentially unitl infinity.
 
 
A more appropriate model would be:
 
 
NoCrash <- function (time, N, parms) {
            with(as.list(parms),
                 {
                   dn <-  r * N * (1 - (N/K)) - h*(N>0)
                   if (! is.numeric(dn)) dn <- 0
                   list(dn)
                 }
                 )
          }

require(deSolve)
parms  = c(r = 0.1, K = 10000, h=200)
init   = c(N = 2000)

out <- ode(func=NoCrash, y = init, parms=parms, times=0:50, method="lsoda")
plot(out)
 
Due to the (N>0) there is only harvesting as long as N>0.
 
note that the plot(out) only works from deSolve version 1.5
 
 
 
Karline Soetaert

________________________________

From: r-sig-dynamic-models-bounces at r-project.org on behalf of Rainer M Krug
Sent: Wed 10/21/2009 12:14 PM
To: r-sig-dynamic-models at r-project.org
Subject: [R-sig-dyn-mod] Problem with simecol and ode - crash from R



Hi

I have a problem with simecol, particularly with the logistic growth model
with harvesting: with the parameters as given below, it crashes R.

When I change h to values below 100, it works, changing the solver to
"euler" or "rk4", it works

Any idea what is happening, or am I doing something wrong?


Thanks,

Rainer

library(simecol)
lgH <- new(
          "odeModel",
          main = function (time, init, parms) {
            with(
                 as.list(c(init, parms)),
                 {
                   dn <-  r * N * (1 - (N/K)) - h
                   list(dn)
                 }
                 )
          },
          parms  = c(r = 0.1, K = 10000, h=200),
          times  = c(from = 0, to = 100, by = 0.5),
          init   = c(N = 2000),
          solver = "lsoda"
)
lgH <-  sim(lgH)
*** stack smashing detected ***: /usr/lib/R/bin/exec/R terminated
======= Backtrace: =========
/lib/tls/i686/cmov/libc.so.6(__fortify_fail+0x48)[0xb7c69138]
/lib/tls/i686/cmov/libc.so.6(__fortify_fail+0x0)[0xb7c690f0]
/home/rkrug/R/i486-pc-linux-gnu-library/2.9/deSolve/libs/deSolve.so[0xb74822b4]
/home/rkrug/R/i486-pc-linux-gnu-library/2.9/deSolve/libs/deSolve.so(dlsoda_+0x1279)[0xb746c329]
/home/rkrug/R/i486-pc-linux-gnu-library/2.9/deSolve/libs/deSolve.so(call_lsoda+0xa2d)[0xb743628d]
/usr/lib/R/lib/libR.so[0xb7d6999c]
/usr/lib/R/lib/libR.so(Rf_eval+0x714)[0xb7d92fb4]
/usr/lib/R/lib/libR.so[0xb7d962fe]
/usr/lib/R/lib/libR.so(Rf_eval+0x451)[0xb7d92cf1]
/usr/lib/R/lib/libR.so[0xb7d94530]
/usr/lib/R/lib/libR.so(Rf_eval+0x451)[0xb7d92cf1]
/usr/lib/R/lib/libR.so(Rf_applyClosure+0x2ac)[0xb7d96c6c]
/usr/lib/R/lib/libR.so(Rf_eval+0x349)[0xb7d92be9]
/usr/lib/R/lib/libR.so[0xb7d23ce6]
/usr/lib/R/lib/libR.so[0xb7dd2d7a]
/usr/lib/R/lib/libR.so(Rf_eval+0x451)[0xb7d92cf1]
/usr/lib/R/lib/libR.so[0xb7d94530]
/usr/lib/R/lib/libR.so(Rf_eval+0x451)[0xb7d92cf1]
/usr/lib/R/lib/libR.so(Rf_applyClosure+0x2ac)[0xb7d96c6c]
/usr/lib/R/lib/libR.so(Rf_eval+0x349)[0xb7d92be9]
/usr/lib/R/lib/libR.so[0xb7d962fe]
/usr/lib/R/lib/libR.so(Rf_eval+0x451)[0xb7d92cf1]
/usr/lib/R/lib/libR.so[0xb7d94530]
/usr/lib/R/lib/libR.so(Rf_eval+0x451)[0xb7d92cf1]
/usr/lib/R/lib/libR.so[0xb7d96487]
/usr/lib/R/lib/libR.so(R_execMethod+0x2c1)[0xb7d96901]
/usr/lib/R/library/methods/libs/methods.so[0xb764b66a]
/usr/lib/R/lib/libR.so[0xb7dd47ab]
/usr/lib/R/lib/libR.so(Rf_eval+0x598)[0xb7d92e38]
/usr/lib/R/lib/libR.so(Rf_applyClosure+0x2ac)[0xb7d96c6c]
/usr/lib/R/lib/libR.so(Rf_eval+0x349)[0xb7d92be9]
/usr/lib/R/lib/libR.so[0xb7d962fe]
/usr/lib/R/lib/libR.so(Rf_eval+0x451)[0xb7d92cf1]
/usr/lib/R/lib/libR.so(Rf_ReplIteration+0x1c5)[0xb7dc0685]
/usr/lib/R/lib/libR.so(run_Rmainloop+0x102)[0xb7dc0a62]
/usr/lib/R/lib/libR.so(Rf_mainloop+0x1c)[0xb7dc0abc]
/usr/lib/R/bin/exec/R(main+0x46)[0x8048776]
/lib/tls/i686/cmov/libc.so.6(__libc_start_main+0xe0)[0xb7b92450]
/usr/lib/R/bin/exec/R[0x8048691]
======= Memory map: ========
08048000-08049000 r-xp 00000000 08:02 811281     /usr/lib/R/bin/exec/R
08049000-0804a000 rw-p 00000000 08:02 811281     /usr/lib/R/bin/exec/R
0804a000-08c8b000 rw-p 0804a000 00:00 0          [heap]
b742f000-b7485000 r-xp 00000000 08:03 10371078
/home/rkrug/R/i486-pc-linux-gnu-library/2.9/deSolve/libs/deSolve.so
b7485000-b7486000 rw-p 00055000 08:03 10371078
/home/rkrug/R/i486-pc-linux-gnu-library/2.9/deSolve/libs/deSolve.so
b7486000-b7487000 rw-p b7486000 00:00 0
b7487000-b74e9000 r-xp 00000000 08:02 1444606
/usr/lib/R/library/stats/libs/stats.so
b74e9000-b74eb000 rw-p 00062000 08:02 1444606
/usr/lib/R/library/stats/libs/stats.so
b74eb000-b750c000 r-xp 00000000 08:02 1810475
/usr/lib/R/library/grDevices/libs/grDevices.so
b750c000-b750d000 rw-p 00020000 08:02 1810475
/usr/lib/R/library/grDevices/libs/grDevices.so
b7555000-b755f000 r-xp 00000000 08:02 1146976    /lib/libgcc_s.so.1
b755f000-b7560000 rw-p 0000a000 08:02 1146976    /lib/libgcc_s.so.1
b7560000-b7562000 r-xp 00000000 08:03 4572127
/home/rkrug/R/i486-pc-linux-gnu-library/2.9/simecol/libs/simecol.so
b7562000-b7563000 rw-p 00001000 08:03 4572127
/home/rkrug/R/i486-pc-linux-gnu-library/2.9/simecol/libs/simecol.so
b7563000-b75d5000 rw-p b7563000 00:00 0
b75d5000-b75d6000 r-xp 00000000 08:02 1226930    /usr/lib/gconv/ISO8859-1.so
b75d6000-b75d8000 rw-p 00000000 08:02 1226930    /usr/lib/gconv/ISO8859-1.so
b75d8000-b7648000 rw-p b75d8000 00:00 0
b7648000-b764e000 r-xp 00000000 08:02 1794882
/usr/lib/R/library/methods/libs/methods.so
b764e000-b764f000 rw-p 00005000 08:02 1794882
/usr/lib/R/library/methods/libs/methods.so
b764f000-b767b000 rw-p b764f000 00:00 0
b767b000-b7684000 r-xp 00000000 08:02 467611     /lib/tls/i686/cmov/
libnss_files-2.7.so
b7684000-b7686000 rw-p 00008000 08:02 467611     /lib/tls/i686/cmov/
libnss_files-2.7.so
b7686000-b768e000 r-xp 00000000 08:02 467613     /lib/tls/i686/cmov/
libnss_nis-2.7.so
b768e000-b7690000 rw-p 00007000 08:02 467613     /lib/tls/i686/cmov/
libnss_nis-2.7.so
b7690000-b76a4000 r-xp 00000000 08:02 467608     /lib/tls/i686/cmov/li
Process R aborted at Wed Oct 21 12:05:59 2009

--
Rainer M. Krug, PhD (Conservation Ecology, SUN), MSc (Conservation Biology,
UCT), Dipl. Phys. (Germany)

Centre of Excellence for Invasion Biology
Natural Sciences Building
Office Suite 2039
Stellenbosch University
Main Campus, Merriman Avenue
Stellenbosch
South Africa

Cell:           +27 - (0)83 9479 042
Fax:            +27 - (0)86 516 2782
Fax:            +49 - (0)721 151 334 888
email:          Rainer at krugs.de

Skype:          RMkrug
Google:         R.M.Krug at gmail.com

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-dynamic-models mailing list
R-sig-dynamic-models at r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-dynamic-models


-------------- next part --------------
A non-text attachment was scrubbed...
Name: not available
Type: application/ms-tnef
Size: 11212 bytes
Desc: not available
URL: <https://stat.ethz.ch/pipermail/r-sig-dynamic-models/attachments/20091021/9face222/attachment.bin>


More information about the R-sig-dynamic-models mailing list