aboutsummaryrefslogtreecommitdiff
path: root/oracles.c
diff options
context:
space:
mode:
authorArun Isaac2021-02-03 13:02:23 +0530
committerArun Isaac2021-02-03 13:02:23 +0530
commit1c9128be78a68805ab0d80631c1984fd4291d894 (patch)
treee0051e36633f502baf559f7f0030aa97ce1fce8c /oracles.c
downloadnsmc-1c9128be78a68805ab0d80631c1984fd4291d894.tar.gz
nsmc-1c9128be78a68805ab0d80631c1984fd4291d894.tar.lz
nsmc-1c9128be78a68805ab0d80631c1984fd4291d894.zip
Initial commit
Diffstat (limited to 'oracles.c')
-rw-r--r--oracles.c141
1 files changed, 141 insertions, 0 deletions
diff --git a/oracles.c b/oracles.c
new file mode 100644
index 0000000..80f401f
--- /dev/null
+++ b/oracles.c
@@ -0,0 +1,141 @@
+#include <math.h>
+#include <gsl/gsl_math.h>
+#include <gsl/gsl_randist.h>
+#include "oracles.h"
+#include "utils.h"
+
+double bernoulli_extent_generator
+(const gsl_rng* r, double p, double r0, double r1)
+{
+ return gsl_ran_bernoulli(r, p) ? r1 : r0;
+}
+
+double bernoulli_true_volume (double p, double r0, double r1, unsigned int dimension)
+{
+ return volume_of_ball(dimension)
+ * (p * gsl_pow_uint(r1, dimension)
+ + (1 - p) * gsl_pow_uint(r0, dimension));
+}
+
+double uniform_extent_generator (const gsl_rng* r, double a, double b)
+{
+ return gsl_ran_flat(r, a, b);
+}
+
+// TODO: Verify the accuracy of this function for non-trivial a, b
+double uniform_true_volume (double a, double b, unsigned int dimension)
+{
+ return exp(ln_volume_of_ball(dimension) + dimension*log(b) - log(dimension + 1))
+ - exp(ln_volume_of_ball(dimension) + dimension*log(a) - log(dimension + 1));
+}
+
+double beta_extent_generator
+(const gsl_rng* r, double alpha, double beta)
+{
+ return gsl_ran_beta(r, alpha, beta);
+}
+
+double beta_true_volume (double alpha, double beta, unsigned int dimension)
+{
+ double vol = volume_of_ball(dimension);
+ for (unsigned int r=0; r<dimension; r++)
+ vol *= (alpha + r) / (alpha + beta + r);
+ return vol;
+}
+
+double symmetric_spiral_extent_oracle (const gsl_vector* x)
+{
+ double angle = atan2(gsl_vector_get(x, 1), gsl_vector_get(x, 0));
+ return fabs(angle);
+}
+
+double right_triangle_extent_oracle (const gsl_vector* x, double base, double height)
+{
+ double max_angle = atan(height / base);
+ double angle = atan2(gsl_vector_get(x, 1), gsl_vector_get(x, 0));
+ return (angle > 0) && (angle < max_angle) ? base / cos(angle) : 0;
+}
+
+double right_triangle_true_volume (double base, double height)
+{
+ return 0.5 * base * height;
+}
+
+double sphere_extent_oracle (const gsl_vector* x, double radius)
+{
+ return radius;
+}
+
+double sphere_maximum_extent (double radius)
+{
+ return radius;
+}
+
+double plane_extent_oracle (const gsl_vector* x, double displacement)
+{
+ return displacement / fabs(gsl_vector_get(x, 0));
+}
+
+double infinity_norm (const gsl_vector* x)
+{
+ double max = fabs(gsl_vector_get(x, 0));
+ for (int i=1; i<x->size; i++)
+ max = GSL_MAX(max, fabs(gsl_vector_get(x, i)));
+ return max;
+}
+
+double cube_extent_oracle (const gsl_vector* x, double edge)
+{
+ return edge / 2 / infinity_norm(x);
+}
+
+double cube_extent_oracle_with_center
+(const gsl_vector* x, const gsl_vector* center, double edge)
+{
+ double min = (0.5*edge - GSL_SIGN(gsl_vector_get(x, 0))*gsl_vector_get(center, 0))
+ / fabs(gsl_vector_get(x, 0));
+ for (int i=1; i<center->size; i++)
+ min = GSL_MIN(min, (0.5*edge - GSL_SIGN(gsl_vector_get(x, i))*gsl_vector_get(center, i))
+ / fabs(gsl_vector_get(x, i)));
+ return min;
+}
+
+double cube_true_volume (double edge, unsigned int dimension)
+{
+ return gsl_pow_uint(edge, dimension);
+}
+
+double cube_maximum_extent (double edge, unsigned int dimension)
+{
+ return edge * sqrt(dimension) / 2;
+}
+
+double ellipsoid_extent_oracle (const gsl_vector* x, const gsl_vector* axes)
+{
+ double k = 0;
+ for (int i=0; i<axes->size; i++)
+ k += gsl_pow_2(gsl_vector_get(x, i) / gsl_vector_get(axes, i));
+ return 1/sqrt(k);
+}
+
+double ellipsoid_true_volume (const gsl_vector* axes)
+{
+ unsigned int dimension = axes->size;
+ double vol = volume_of_ball(dimension);
+ for (int i=0; i<dimension; i++)
+ vol *= gsl_vector_get(axes, i);
+ return vol;
+}
+
+double spheroid_extent_oracle (const gsl_vector* x, double eccentricity)
+{
+ unsigned int dimension = x->size;
+ gsl_vector_const_view xsub = gsl_vector_const_subvector(x, 1, dimension - 1);
+ return 1/sqrt(gsl_pow_2(gsl_blas_dnrm2(&xsub.vector))
+ + gsl_pow_2(gsl_vector_get(x, 0) / eccentricity));
+}
+
+double spheroid_true_volume (double eccentricity, unsigned int dimension)
+{
+ return volume_of_ball(dimension) * eccentricity;
+}