Commit 4682d721 authored by Shawn Samuel Carl McAdam's avatar Shawn Samuel Carl McAdam
Browse files

comments and more

parent 3575a541
......@@ -49,6 +49,7 @@ int main(int argc, char *argv[]){
int char_walker=1;
char current_char=current_arg[char_walker];
while(current_char != '\0'){
// TODO switch case
if(current_char == 'h')
usage(EXIT_SUCCESS, "");
else if(current_char == 'a'){
......
#include <inttypes.h> // uint64_t
#include "bin_dfs.h"
// private helper functions
long double cdf_helper(uint64_t k, uint64_t n, long double p, long double sum, long double error);
long double one_minus_cdf_helper(uint64_t k, uint64_t n, long double p, long double sum, long double error);
long double pmf_helper(uint64_t k, uint64_t n, long double p);
// declarations of public driver functions
long double pmf(uint64_t k, uint64_t n, long double p){
// TODO store even numbers as private args to helper functions and bitshift at the end
// probably could store a pointer to our even numbers...?
long double bin_pmf(uint64_t k, uint64_t n, long double p){
if(k > n)
return 0.0; // the function is being used improperly so return 0.0.
if(2*k>n && p > 0.5)
......@@ -16,44 +11,43 @@ long double pmf(uint64_t k, uint64_t n, long double p){
return pmf_helper(k, n, p);
}
long double cdf(uint64_t k, uint64_t n, long double p){
if(2*k>n)
return 1-one_minus_cdf_helper(k+1, n, p, 0.0, 0.0);
return cdf_helper(k, n, p, 0.0, 0.0);
}
// Can't do tail recursion here because the intermediate products can get
// far larger than LDBL_MAX
long double pmf_helper(uint64_t k, uint64_t n, long double p){
long double bin_pmf_helper(uint64_t k, uint64_t n, long double p){
// Can't do tail recursion here because the intermediate products can get
// far larger than LDBL_MAX
// work through our k's then work through our n's.
if(n == 0)
return 1.0;
if(k == 0)
return (1-p)*pmf_helper(0, n-1, p);
return n*p*pmf_helper(k-1, n-1, p)/k;
return n*p*bin_pmf_helper(k-1, n-1, p)/k;
}
long double cdf(uint64_t k, uint64_t n, long double p){
if(2*k>n)
return 1-one_minus_cdf_helper(k+1, n, p, 0.0, 0.0);
return cdf_helper(k, n, p, 0.0, 0.0);
}
// accumulate a sum and error using tail recursion. Call with sum = error = 0.0
long double cdf_helper(uint64_t k, uint64_t n, long double p, long double sum, long double error){
long double bin_cdf_helper(uint64_t k, uint64_t n, long double p, long double sum, long double error){
// sum until k == 0
if(k == 0)
return sum + pmf(0, n, p) + error; // add error back to sum to get
return sum + pmf(0, n, p) + error; // add error back to sum to get a more accurate answer
// perform a compensated addition on sum
long double temp = pmf(k, n, p);
long double next = sum + temp;
return cdf_helper(k-1, n, p, next, error+((sum-next)+temp));
return bin_cdf_helper(k-1, n, p, next, error+((sum-next)+temp));
}
// calculates 1-cdf(k). Can result in fewer steps to find a solution to 1-cdf(k) if 2*k>n.
long double one_minus_cdf_helper(uint64_t k, uint64_t n, long double p, long double sum, long double error){
long double bin_one_minus_cdf_helper(uint64_t k, uint64_t n, long double p, long double sum, long double error){
// sum until k == n
if(k > n)
return sum + error; // add error back to sum to get
return sum + error; // add error back to sum to get a more accurate answer
// perform a compensated addition on sum
long double temp = pmf(k, n, p);
long double next = sum + temp;
return one_minus_cdf_helper(k+1, n, p, next, error+((sum-next)+temp));
return bin_one_minus_cdf_helper(k+1, n, p, next, error+((sum-next)+temp));
}
......@@ -12,7 +12,7 @@
* Returns:
* (n_choose_k)*(p^k)*((1-p)^(n-k)) if k<=n, 0 otherwise
*/
long double pmf(uint64_t k, uint64_t n, long double p);
long double bin_pmf(uint64_t k, uint64_t n, long double p);
/*
......@@ -24,7 +24,12 @@ long double pmf(uint64_t k, uint64_t n, long double p);
* P(X<=k) where X~bin(n, p)
* ie, \sum_{j=0}^k pmf(j, n, p)
*/
long double cdf(uint64_t k, uint64_t n, long double p);
long double bin_cdf(uint64_t k, uint64_t n, long double p);
// private helper functions called by pmf and cdf
long double bin_cdf_helper(uint64_t k, uint64_t n, long double p, long double sum, long double error);
long double bin_one_minus_cdf_helper(uint64_t k, uint64_t n, long double p, long double sum, long double error);
long double bin_pmf_helper(uint64_t k, uint64_t n, long double p);
#endif // DIST_FUNCT_H
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment