R Under development (unstable) (2024-05-20 r86569 ucrt) -- "Unsuffered Consequences"
Copyright (C) 2024 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> if (!requireNamespace("spatstat.geom"))
+     q("no")
Loading required namespace: spatstat.geom
> library("polyCub")
> 
> ## bivariate, isotropic Gaussian density
> f <- function (s, mean, sd)
+     dnorm(s[,1], mean=mean[1], sd=sd) * dnorm(s[,2], mean=mean[2], sd=sd)
> 
> ## circular domain represented by a polygon
> r <- 5
> center <- c(3,2)
> npoly <- 128
> disc.owin <- spatstat.geom::disc(radius=r, centre=center, npoly=npoly)
> 
> ## parameters for f
> m <- c(1,1)
> sd <- 3
> 
> ## target value of the integral over the _polygonal_ circle
> intExact <- 0.65844436  # taken from exact.Gauss cubature (below)
> stopIfDiff <- function(int, ...)
+     if(!isTRUE(all.equal.numeric(intExact, int, ..., check.attributes = FALSE))) {
+         if (is.call(cl <- substitute(int))) cl <- cl[1]
+         stop(deparse(cl), " result not equal to reference value")
+     }
> 
> ## reproduce saved reference value
> if (identical(Sys.getenv("R_GPCLIBPERMIT"), "true") &&
+     local({pkg <- "gpclib"; requireNamespace(pkg)}) && # undeclared ...
+     requireNamespace("mvtnorm"))
+     stopIfDiff(polyCub.exact.Gauss(disc.owin, mean=m, Sigma=sd^2*diag(2)),
+                tolerance = 1e-8)
> 
> ## exact value of the integral over the _real_ circle
> stopIfDiff(circleCub.Gauss(center=center, r=r, mean=m, sd=sd),
+            tolerance = 0.001)  # agreement depends on 'npoly'
> 
> ## polyCub.midpoint
> stopIfDiff(polyCub.midpoint(disc.owin, f, mean=m, sd=sd, dimyx=500),
+            tolerance = 0.001)
> 
> ## polyCub.SV
> intC <- polyCub.SV(disc.owin, f, mean=m, sd=sd, nGQ=3, engine="C")
> intR <- polyCub.SV(disc.owin, f, mean=m, sd=sd, nGQ=3, engine="R")
> stopifnot(all.equal(intC, intR))
> stopIfDiff(intC, tolerance = 0.0001)
> 
> ## polyCub.iso (using a numerical approximation of intrfr)
> stopIfDiff(polyCub.iso(disc.owin, f, mean=m, sd=sd, center=m))
> 
> proc.time()
   user  system elapsed 
   1.23    0.06    1.28