In an introductory probability class, one computes the probability of
getting all of n possible coupons in r individual purchases. The naive
approach with inclusion-exclusion leads to the awful formula
f(n,r) = \sum_{k=0}^n \binom{n,k} \frac{(n-k)^r}{n^r}
Computing this in GP is straightforward:
(10:32) gp > g(n,k,r)=(-1)^(k)*binomial(n,k)*(n-k)^r/n^r
%3 = (n,k,r)->(-1)^(k)*binomial(n,k)*(n-k)^r/n^r
(16:43) gp > f(n,r)=sum(k=0,n,1.0*g(n,k,r))
%4 = (n,r)->sum(k=0,n,1.0*g(n,k,r))
(16:43) gp > f(6,12)
%5 = 0.43781568062117902081322291656082236787
(16:43) gp > f(50,100)
%6 = 0.00016616318861823418331310572392871711604
(16:45) gp > f(50,200)
%7 = 0.39828703196689437232172533821416865689
(16:47) gp > f(365,2236)
Doing it in SageMath seems to be much more delicate. This, for example,
fails:
sage: g(n,k,r)=(-1)^(k)*binomial(n,k)*(n-k)^r/n^r
sage: f(n,r)=sum(1.0*g(n,k,r),k,0,n)
sage: f(365,2000).unhold()
0.0
It works better if I rewrite the function with ((n-k)/n)^r, of course. I
presume the issue is that the default precision for real numbers is not
enough to handle the computation. How do I tell SageMath to use a higher
precision?
Fernando
--
=============================================================
Fernando Q. Gouveahttp://www.colby.edu/~fqgouvea
Carter Professor of Mathematics
Dept. of Mathematics
Colby College
5836 Mayflower Hill
Waterville, ME 04901
Reached by phone, a Westbury woman who identified herself as Onorato's
granddaughter said he was a retired businessman. He lived alone and
never married or had children, she said.
-- New York Newsday, December 27, 2003
--
You received this message because you are subscribed to the Google Groups
"sage-support" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to sage-support+unsubscr...@googlegroups.com.
To view this discussion on the web visit
https://groups.google.com/d/msgid/sage-support/fd3261d8-3ee5-48a8-ba6e-40530790b92f%40colby.edu.
