[Rd] uniroot violates bounds?

[Serguei Sokol]

Le 18/02/2023 à 21:44, J C Nash a écrit :
> I wrote first cut at unirootR for Martin M and he revised and put in
> Rmpfr.
>
> The following extends Ben's example, but adds the unirootR with trace
> output.
>
> c1 <- 4469.822
> c2 <- 572.3413
> f <- function(x) { c1/x - c2/(1-x) }; uniroot(f, c(1e-6, 1))
> uniroot(f, c(1e-6, 1))
> library(Rmpfr)
> unirootR(f, c(1e-6, 1), extendInt="no", trace=1)
>
> This gives more detail on the iterations, and it looks like the Inf is 
> the
> issue. But certainly we could do more to avoid "gotchas" like this. If
> someone is interested in some back and forth, I'd be happy to give it a
> try, but I think progress would be better with more than one contributor.
For me, the following fix makes the job :

--- nlm.R.old    2018-09-25 10:44:49.000000000 +0200
+++ nlm.R    2023-02-20 10:46:39.893542531 +0100
@@ -143,14 +143,14 @@

      if(check.conv) {
      val <- tryCatch(.External2(C_zeroin2, function(arg) f(arg, ...),
-                   lower, upper, f.lower, f.upper,
+                   lower, upper, f.low., f.upp.,
                     tol, as.integer(maxiter)),
              warning = function(w)w)
      if(inherits(val, "warning"))
          stop("convergence problem in zero finding: ", 
conditionMessage(val))
      } else {
      val <- .External2(C_zeroin2, function(arg) f(arg, ...),
-              lower, upper, f.lower, f.upper,
+              lower, upper, f.low., f.upp.,
                tol, as.integer(maxiter))
      }
      iter <- as.integer(val[2L])


Best,
Serguei.

>
> Best,
>
> John Nash
>
> On 2023-02-18 12:28, Ben Bolker wrote:
>> c1 <- 4469.822
>> c2 <- 572.3413
>> f <- function(x) { c1/x - c2/(1-x) }; uniroot(f, c(1e-6, 1))
>> uniroot(f, c(1e-6, 1))
>>
>>
>>     provides a root at -6.00e-05, which is outside of the specified 
>> bounds.  The default value of the "extendInt" argument to uniroot() 
>> is "no", as far as I can see ...
>>
>> $root
>> [1] -6.003516e-05
>>
>> $f.root
>> [1] -74453981
>>
>> $iter
>> [1] 1
>>
>> $init.it
>> [1] NA
>>
>> $estim.prec
>> [1] 6.103516e-05
>>
>>
>>    I suspect this fails because f(1) (value at the upper bound) is 
>> infinite, although setting interval to c(0.01, 1) does work/give a 
>> sensible answer ...  (works for a lower bound of 1e-4, fails for 1e-5 
>> ...)
>>
>>    Setting the upper bound < 1 appears to avoid the problem.
>>
>>   For what it's worth, the result has an "init.it" component, but the 
>> only thing the documentation says about it is " component ‘init.it’ 
>> was added in R 3.1.0".
>>
>>    And, I think (?) that the 'trace' argument only produces any 
>> output if the 'extendInt' option is enabled?
>>
>>    Inspired by 
>> https://stackoverflow.com/questions/75494696/solving-a-system-of-non-linear-equations-with-only-one-unknown/75494955#75494955
>>
>>    cheers
>>     Ben Bolker
>>
>> ______________________________________________
>> R-devel using r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-devel
>
> ______________________________________________
> R-devel using r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel