I am running uniroot
function in the loop, but hit an error and the code stopped. The code is as below;
func <-function(f) -a*b/c*0.5*d*e^2 + (d/f-1)*g*sin(h*(pi/180))-i
dat <- data.frame(a = c(0.99,0.99,0.99),
b = c(0.1986572,0.1986572,0.1986572),
c = c(237.5,237.5,237.5),
d = c(1028.372, 1028.711, 1028.372),
e = c(2.46261, 2.986461, 2.46261),
f = c(-1,-1,-1),
g = c(9.8,9.8,9.8),
h = c(-54.97964, -51.65978, -54.97964),
i = c(0.03699588, -0.0375189, 0.03699588))
for(j in 1:length(dat$a)){
a <- dat$a[j]
b <- dat$b[j]
c <- dat$c[j]
d <- dat$d[j]
e <- dat$e[j]
#f: this should be solved by uniroot
g <- dat$g[j]
h <- dat$h[j]
i <- dat$i[j]
sol <- uniroot(func,c(0, 2000),extendInt = "yes")
dat$f[j] <- sol$root
print(j)
}
Running above code, hit the below error:
[1] 1
Error in uniroot(func, c(0, 2000), extendInt = "yes") :
no sign change found in 1000 iterations
The code stopped at j=1
, and did not go to j=2 & 3
. Therefore, dat$f
shows
> dat$f
[1] 1526.566 -1.000 -1.000
My goal is when uniroot
hits an error in a given j
, put NA
in dat$f[j]
, and continue the loop by the end.
If this works, dat$f[1]
and dat$f[3]
should have the same value (=1526.566) using the the above dataframe.
Please advise me on how to deal with the uniroot
error.
The code in the question would work if the lower bound is set to 1 instead of 0. The problem is that if f is 0 then func is undefined since f is in the denominator.
While that is sufficient it is suggested you make these changes:
try
to allow uniroot to continue to run even if there are errors (although there are none in t his example)Use the more compact form of func shown and put the results in res.
func <- function(f, data) with(data, -a*b/c*0.5*d*e^2 + (d/f-1)*g*sin(h*(pi/180))-i)
nr <- nrow(dat)
res <- numeric(nr)
for(j in 1:nr){
sol <- try(uniroot(func, c(1, 2000), data = dat[j, ], extendInt = "yes") )
res[j] <- if (inherits(sol, "try-error")) NA else sol$root
print(j)
}
## [1] 1
## [1] 2
## [1] 3
res
## [1] 1526.566 2014.476 1526.566
dat <- data.frame(a = c(0.99,0.99,0.99),
b = c(0.1986572,0.1986572,0.1986572),
c = c(237.5,237.5,237.5),
d = c(1028.372, 1028.711, 1028.372),
e = c(2.46261, 2.986461, 2.46261),
g = c(9.8,9.8,9.8),
h = c(-54.97964, -51.65978, -54.97964),
i = c(0.03699588, -0.0375189, 0.03699588))