numer-utils: Numerical Utilities - Functions, Constants

Description

The DPQ package provides some numeric constants used in some of its distribution computations.

all_mpfr() and any_mpfr() return TRUE iff all (or ‘any’, respectively) of their arguments inherit from class "mpfr" (from package Rmpfr).

logr(x,a) computes log(x / (x + a)) in a numerically stable way.

Usage

## Numeric Constants : % mostly in   ../R/beta-fns.R
M_LN2        # = log(2)  = 0.693....
M_SQRT2      # = sqrt(2) = 1.4142...
M_cutoff     # := If |x| > |k| * M_cutoff, then  log[ exp(-x) * k^x ]  =~=  -x
             #  = 3196577161300663808 ~= 3.2e+18
M_minExp     # = log(2) * .Machine$double.min.exp # ~= -708.396..
G_half       # = sqrt(pi) = Gamma( 1/2 )
## Functions :
all_mpfr(...)
any_mpfr(...)
logr(x, a)    # == log(x / (x + a)) -- but numerically smart; x >= 0, a > -x
okLongDouble(lambda = 999, verbose = 0L, tol = 1e-15)

Arguments

…

numeric or "mpfr" numeric vectors.

x, a

number-like, not negative, now may be vectors of length(.) > 1.

lambda

a number, typically in the order of 500--10'000.

verbose

a non-negative integer, if not zero, okLongDouble() prints the intermediate long double computations' results.

tol

numerical tolerance used to determine the accuracy required for near equality in okLongDouble().

Value

The numeric constant in the first case; a numeric (or "mpfr") vector of appropriate size in the 2nd case.

okLongDouble() returns a logical, TRUE iff the long double arithmetic with expl() and logl() seems to work accurately

Details

all_mpfr(),
all_mpfr() :: test if all or any of their arguments or of class "mpfr" (from package Rmpfr). The arguments are evaluated only until the result is determined, see the example.
logr(): computes \(\log( x / (x+a) )\) in a numerically stable way.

Examples

Run this code

# NOT RUN {
(Ms <- ls("package:DPQ", pattern = "^M"))
lapply(Ms, function(nm) { cat(nm,": "); print(get(nm)) }) -> .tmp

logr(1:3, a=1e-10)

okLongDouble() # typically TRUE, but not e.g. in a valgrinded R-devel of Oct.2019
## Here is typically the "boundary":
okLongDouble(11355, verbose=TRUE) # typically TRUE  (also for lambda <= 11355)
okLongDouble(11356, verbose=TRUE) # typically FALSE (also for lambda >= 11356)
# }