[R-sig-dyn-mod] bvp chemical reactor?

Radovan Omorjan omorr at uns.ac.rs
Fri May 6 13:56:38 CEST 2016


Also, guessing values will also help as well (not using the absolute 
value of y[1]). I just gave here these rather poor guess values:
-first state variable given 1, and second state variable given 0

xguess <- seq(0,1,len=50)
yguess <- rbind(rep(1,50),rep(0,50))

and both bvptwp() and bvpcol() will work correctly giving the expected 
result .

Regards,
Radovan

On 5/5/2016 3:03 PM, Radovan Omorjan wrote:
> I found out that if I define the function this way (use absolute value 
> of y[1]) then the result will be correct and both solvers will work.
> Could someone try to explain this to me please.
>
> Best Regards,
> Radovan
>
> THIS WILL WORK!
> -----------------------
> fun <- function(x, y, p) {
>   dy1 <- y[2]
>   dy2 <- Pe*Da*abs(y[1])^n + Pe*y[2]
>   return(list(c(dy1,dy2)))
> }
> --------------------------
> On 5/5/2016 2:51 PM, Radovan Omorjan wrote:
>> Hello,
>>
>> This one is also one of my standard examples (chemical reactor) of 
>> two point boundary problem.
>> I tried the three functions (bvptwp(), bvpcol(), bvpshoot()) with the 
>> default parameters and only bvpshoot() gave the results.
>> Many other examples gave me no problem in solving them with all the 
>> three solvers, but this particular one.
>>
>> QUESTION: I would be very grateful if someone could take a look at 
>> this example and give me a hint if there is a way to get the result 
>> with the other two solvers.
>> Just to mention that I do not have to much experience with these 
>> solvers and I do not know how to explain my students why these two do 
>> not work for this particular example.
>>
>> Best Regards,
>> Radovan
>>
>> ---------------------------
>> require(bvpSolve)
>>
>> fun <- function(x, y, p) {
>>   dy1 <- y[2]
>>   dy2 <- Pe*Da*y[1]^n + Pe*y[2]
>>   return(list(c(dy1,dy2)))
>> }
>>
>> Pe <- 1
>> Da <- 2
>> n <- 0.5
>>
>> # initial and ending conditions
>> init <- c(y <- 1, dy <- NA)
>> end <- c(y <- NA, dy <- 0)
>>
>> # Solve bvp
>> #
>> # This will not work
>> #
>> # sol  <- bvptwp(yini = init, x = seq(0,1,len=50),
>> #               func = fun, yend = end)
>> # Error in bvpsolver(1, yini, x, func, yend, parms, order, ynames, 
>> xguess,  :
>> #                     The Expected No. Of mesh points Exceeds Storage 
>> Specifications.
>> # sol  <- bvpcol(yini = init, x = seq(0,1,len=50),
>> #               func = fun, yend = end)
>> # Error in bvpsolver(2, yini, x, func, yend, parms, order, ynames, 
>> xguess,  :
>> #                     The collocation matrix is singular for the 
>> final continuation problem
>> # This one will work
>> #
>> sol  <- bvpshoot(yini = init, x = seq(0,1,len=50),
>>                func = fun, yend = end, guess=-1)
>>
>> x <- sol[,1]
>> y <- sol[,2]
>> par(mar=c(4,4,1.5,1.5),mex=.8,mgp=c(2,.5,0),tcl=0.3)
>> plot(x,y)
>>
>> _______________________________________________
>> R-sig-dynamic-models mailing list
>> R-sig-dynamic-models at r-project.org
>> https://stat.ethz.ch/mailman/listinfo/r-sig-dynamic-models
>>
>
> _______________________________________________
> R-sig-dynamic-models mailing list
> R-sig-dynamic-models at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-dynamic-models
>



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