1 //===----------------------------------------------------------------------===//
2 //
3 // Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
4 // See https://llvm.org/LICENSE.txt for license information.
5 // SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
6 //
7 //===----------------------------------------------------------------------===//
8 //
9 // REQUIRES: long_tests
10
11 // <random>
12
13 // template<class RealType = double>
14 // class fisher_f_distribution
15
16 // template<class _URNG> result_type operator()(_URNG& g);
17
18 #include <random>
19 #include <cassert>
20 #include <vector>
21 #include <algorithm>
22 #include <cmath>
23
24 #include "test_macros.h"
25
fac(double x)26 double fac(double x)
27 {
28 double r = 1;
29 for (; x > 1; --x)
30 r *= x;
31 return r;
32 }
33
34 double
I(double x,unsigned a,unsigned b)35 I(double x, unsigned a, unsigned b)
36 {
37 double r = 0;
38 for (int j = a; static_cast<unsigned>(j) <= a+b-1; ++j)
39 r += fac(a+b-1)/(fac(j) * fac(a + b - 1 - j)) * std::pow(x, j) *
40 std::pow(1-x, a+b-1-j);
41 return r;
42 }
43
44 double
f(double x,double m,double n)45 f(double x, double m, double n)
46 {
47 return I(m * x / (m*x + n), static_cast<unsigned>(m/2), static_cast<unsigned>(n/2));
48 }
49
main(int,char **)50 int main(int, char**)
51 {
52 // Purposefully only testing even integral values of m and n (for now)
53 {
54 typedef std::fisher_f_distribution<> D;
55 typedef std::mt19937 G;
56 G g;
57 D d(2, 4);
58 const int N = 100000;
59 std::vector<D::result_type> u;
60 for (int i = 0; i < N; ++i)
61 {
62 D::result_type v = d(g);
63 assert(v >= 0);
64 u.push_back(v);
65 }
66 std::sort(u.begin(), u.end());
67 for (int i = 0; i < N; ++i)
68 assert(std::abs(f(u[i], d.m(), d.n()) - double(i)/N) < .01);
69 }
70 {
71 typedef std::fisher_f_distribution<> D;
72 typedef std::mt19937 G;
73 G g;
74 D d(4, 2);
75 const int N = 100000;
76 std::vector<D::result_type> u;
77 for (int i = 0; i < N; ++i)
78 {
79 D::result_type v = d(g);
80 assert(v >= 0);
81 u.push_back(v);
82 }
83 std::sort(u.begin(), u.end());
84 for (int i = 0; i < N; ++i)
85 assert(std::abs(f(u[i], d.m(), d.n()) - double(i)/N) < .01);
86 }
87 {
88 typedef std::fisher_f_distribution<> D;
89 typedef std::mt19937 G;
90 G g;
91 D d(18, 20);
92 const int N = 100000;
93 std::vector<D::result_type> u;
94 for (int i = 0; i < N; ++i)
95 {
96 D::result_type v = d(g);
97 assert(v >= 0);
98 u.push_back(v);
99 }
100 std::sort(u.begin(), u.end());
101 for (int i = 0; i < N; ++i)
102 assert(std::abs(f(u[i], d.m(), d.n()) - double(i)/N) < .01);
103 }
104
105 return 0;
106 }
107