Get best of the week

The use of R in precision calculations

Periodically, there are problems, even in everyday life, when the bit accuracy float64/ int64is not enough to get an answer with the required accuracy. Rushing about in search of another instrument? Also an option.


Or you can not do this, but be curious and find out that for calculating with arbitrary accuracy the GNU MPFR library has long been made to which there are wrappers for almost all languages. Practice shows that few people are familiar with this library, which is probably caused by the peculiarities of study programs at universities and the subsequent programmer mainstream.


The library is good and deserves to be paid attention to, at least within the framework of broadening its horizons. On R to it there is an Rmpfr wrapper . Below I will give a simple example on problems for schoolchildren (well, do not touch the design data under the NDA) and I will touch on a number of classic rakes that are attacked almost immediately.


It is a continuation of previous publications .


For example, take problem No. 39 from the 2019/2020 school year competition in the Quantum: The
positive numbers x and y are such that in the inequality below the left fraction is larger than the right. Which is greater: x or y?


chain fractions


Naturally, it must be solved analytically (alternating inequality), but this is not a reason not to use it as a demonstration.


We draw a simple code to calculate the sequence of fractions. You can immediately set the final value, but then the demoscene for mpfr, so we move in small steps.


We solve by standard methods
library(tidyverse)
library(magrittr)
options(digits=15)

frec <- function(stopval, n, x){
  res <- ifelse(n == stopval, 
                (stopval - 1) + stopval/(stopval + 1 + x), 
                (n - 1 ) + n / (frec(stopval, n + 2, x))
  )
  res
}

frec_wrap <- function(stopval, x){
  res <- frec(stopval = stopval, n = 1, x = x)
  print(glue::glue("{stopval}: {res}"))
  res
}

sol_df <- tibble(stopval = seq(1, 29, by = 2)) %>%
  mutate(val1 = purrr::map_dbl(stopval, frec_wrap, x = 1),
         val2 = purrr::map_dbl(stopval, frec_wrap, x = 5),
         delta = val1 - val2)

And here is a bad luck, already at the 14th iteration (stop number = 29), accuracy is not enough for us to distinguish between fractions. And you have to consider right up to 2019!


Failure 1


? , Rmfpr. — float64 mpfr .


library(tidyverse)
library(magrittr)
library(Rmpfr)
frec2 <- function(stopval, n, x){
  if(n == stopval){
    (stopval - 1) + stopval/(stopval + 1 + x)} else {
      (n - 1 ) + n / (frec2(stopval, n + 2, x))
    }
}
frec2_wrap <- function(stopval, x){
  .precision <- 5000
  res <- frec2(stopval = mpfr(stopval, .precision), 
               n = mpfr(1, .precision), 
               x = mpfr(x, .precision)
  )
  print(glue::glue("{stopval}: {formatMpfr(res, drop0trailing = TRUE)}"))
  res
}

sol2_df <- tibble(stopval = seq(1, 29, by = 2)) %>%
  mutate(val1 = purrr::map(stopval, frec2_wrap, x = 1),
         val2 = purrr::map(stopval, frec2_wrap, x = 5))

. , tibble . .


Failure 2


№1:
tibble , vctrs. ( mpfr S4 ) list-column. - .
, vctrs . , :


for(jj in 1:12){
  #   ,       
  flags_df[[glue("bp_2_{jj}_T")]] <- flags_df[[glue("bp_2_{jj}_in")]] &   flags_df[[glue("flag_2_{jj}")]]
  flags_df[[glue("bp_2_{jj}_F")]] <- flags_df[[glue("bp_2_{jj}_in")]] & ! flags_df[[glue("flag_2_{jj}")]]
}

№2:
list-column . mpfr list-column.


№3
rpfm . StackOverflow . . R?


? --! R, Python! !
. RTFM, .


  • №1. — tibble . data.frame.
  • №2. rpfm . , . data.frame .
  • №3. formatMpfr .

.


library(tidyverse)
library(magrittr)
library(Rmpfr)
#    
frec2 <- function(stopval, n, x){
  if(n == stopval){
    (stopval - 1) + stopval/(stopval + 1 + x)} else {
      (n - 1 ) + n / (frec2(stopval, n + 2, x))
    }
}
frec2_wrap <- function(stopval, x){
  # browser()
  .precision <- 5000
  res <- frec2(stopval = mpfr(stopval, .precision), 
               n = mpfr(1, .precision), 
               x = mpfr(x, .precision)
  )
  print(glue::glue("{stopval}: {formatMpfr(res, drop0trailing = TRUE)}"))
  res
}

sol_df <- data.frame(stopval = seq(111, 119, by = 2)) %>%
  #    mpfr   
  mutate(val1 = new("mpfr", unlist(purrr::map(stopval, frec2_wrap, x = 1))),
         val2 = new("mpfr", unlist(purrr::map(stopval, frec2_wrap, x = 5))),
         delta = val1 - val2)

sol_txt_df <- sol_df %$%
  tibble(stopval = stopval,
         val1_txt = formatMpfr(val1, drop0trailing = TRUE),
         val2_txt = formatMpfr(val2, drop0trailing = TRUE),
         delta_txt = formatMpfr(delta, drop0trailing = TRUE))

P.S. . , . , , GNU MPFR.


— « R ?».

More articles:


All Articles