[R] Fwd: uniroot problem

[Troels Ring]

-------- Videresendt meddelelse --------
Emne: 	Re: [R] uniroot problem
Dato: 	Fri, 16 Sep 2022 10:05:43 +0200
Fra: 	Troels Ring <tring using gvdnet.dk>
Til: 	Andrew Simmons <akwsimmo using gmail.com>



Thanks a lot  - sorry for the confusion - pH[1] was the known pH 4.57.

  But you were rigth: I made simple errors:

ASTA was exactly as you deduced

ASTA <- H + SID - kw/H - ATOT*KA/(KA+H)

And this cures the problem! And had nothing to do with uniroot!
Thanks a lot!!!

Den 16-09-2022 kl. 08:21 skrev Andrew Simmons:
> In this line:
>
> for (i in 1:100) {
>      ATOT[i] <- uniroot(ASTA,c(-1e3,1e3),tol=1e-20,maxiter=1e4,
>                         SID=SID,KA=KA[i],H=10^-pH[1])$root
> }
>
> you've used variable pH, but have no definition for this variable.
>
>
> Also, I can't tell which location has a mistake, but these aren't equivalent:
>
> ASTA <- H + SID - kw/H - ATOT*KA/KA+H
> ATOT <-  (H+SID-kw/H)*(KA+H)/KA
>
> I'm assuming the first one is incorrect, it seems redundant to
> multiply then divide by KA, so change ASTA to:
>
> ASTA <- function(SID,H,ATOT,KA) H + SID - kw/H - ATOT*KA/(KA+H)
>
> I tried something like this:
>
> ASTA <- function(SID,H,ATOT,KA) H + SID - kw/H - ATOT*KA/(KA+H)
>
> SID <- 0.01199
> H <- 10^-4.57
> kw <- 1e-14
>
> KA <- 10^-seq(1, 10, length.out = 100)
>
> ATOT1 <- (H + SID - kw/H) * (KA + H)/KA
>
> ATOT2 <- vapply(KA, function(ka) {
>       uniroot(ASTA, c(-10000, 10000), SID = SID, H = H, KA = ka, tol =
> 1e-20, maxiter = 10000)$root
> }, 0)
>
>
> ASTA(SID, H, ATOT1, KA)
> ASTA(SID, H, ATOT2, KA)
>
> and it seems to work for me. The biggest absolute diff is ~1e-13
>
> On Fri, Sep 16, 2022 at 2:05 AM Troels Ring<tring using gvdnet.dk>  wrote:
>> Dear friends - I have a problem with root finding using uniroot as goes:
>>
>> I have the function
>>
>> ASTA <- function(SID,H,ATOT,KA) {H + SID - kw/H - ATOT*KA/KA+H }
>>
>> I know SID = 0.01199 - and I know H = 10^-4.57 - and I know kw = 1e-14
>>
>> Therefore, as I make KA = 10^-seq(1,10,length=100), by letting ASTA = 0,
>> I can find ATOT for each KA as
>>
>> ATOT <-  (H+SID-kw/H)*(KA+H)/KA
>>
>>
>> and
>>
>>      summary(ATOT) =    Min.  1st Qu.   Median     Mean  3rd Qu. Max.
>>                                         0.012    0.013    0.115 171.263
>> 18.278 3234.407
>>
>> However, if I make use of uniroot by:
>>
>> KA <- seq(1,10,length=100)
>> KA <- 10^-KA
>>
>> ATOT <- c()
>>
>> for (i in 1:100) {
>>      ATOT[i] <- uniroot(ASTA,c(-1e3,1e3),tol=1e-20,maxiter=1e4,
>>                         SID=SID,KA=KA[i],H=10^-pH[1])$root
>> }
>>
>> I get 100 identical ATOT values of 0.01204383 and no warnings.
>>
>> I'm on Windows 10, using slightly old
>>
>> $version.string
>> [1] "R version 4.0.5 (2021-03-31)"
>>
>> Best wishes
>>
>> Troels Ring, MD
>> Aalborg, Denmark
>>
>> ______________________________________________
>> R-help using r-project.org  mailing list -- To UNSUBSCRIBE and more, see
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guidehttp://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
> ______________________________________________
> R-help using r-project.org  mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guidehttp://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
	[[alternative HTML version deleted]]