diff --git a/MailLib/api/com.upokecenter.mail.MediaType.html b/MailLib/api/com.upokecenter.mail.MediaType.html index 13d5ba65a..cc95e9c41 100755 --- a/MailLib/api/com.upokecenter.mail.MediaType.html +++ b/MailLib/api/com.upokecenter.mail.MediaType.html @@ -258,7 +258,7 @@

isMultipart

getParameters

-
public final Map<String,String> getParameters() Gets a list of the parameter names contained in this media type object and  their values. Each parameter name will be in lower case; that is, with its  basic uppercase letters ("A" to "Z") converted to basic lowercase letters  ("a" to "z").
+
public final Map<String,String> getParameters() Gets a list of the parameter names contained in this media type object and  their values. Each parameter name will be in lowercase; that is, with its  basic uppercase letters ("A" to "Z") converted to basic lowercase letters  ("a" to "z").
 

Returns:

diff --git a/MailLib/api/com.upokecenter.text.Idna.html b/MailLib/api/com.upokecenter.text.Idna.html index f63fec004..b445b373c 100755 --- a/MailLib/api/com.upokecenter.text.Idna.html +++ b/MailLib/api/com.upokecenter.text.Idna.html @@ -33,11 +33,10 @@

com.upokecenter.text.Idna

requirements for labels with such characters. An example is "eá".

An A-label is an LDH label beginning with "xn--" where the letters can be any - combination of basic uppercase and/or basic lowercase letters, and is + combination of basic uppercase and basic lowercase letters, and is convertible to a U-label. An example is "xn--e-ufa".

An XN-label is - an LDH label beginning with "xn--" where the letters can be any combination - of basic uppercase and/or basic lowercase letters.

+ an LDH label beginning with "xn--" where the letters can be any combination of basic uppercase and basic lowercase letters.

NOTICE: While this class's source code is in the public domain, the class uses two internal classes, called NormalizationData and IdnaData, that include diff --git a/MailLib/api/com.upokecenter.text.NormalizerInput.html b/MailLib/api/com.upokecenter.text.NormalizerInput.html index f3354ec81..2c39aafd9 100755 --- a/MailLib/api/com.upokecenter.text.NormalizerInput.html +++ b/MailLib/api/com.upokecenter.text.NormalizerInput.html @@ -19,7 +19,7 @@

com.upokecenter.text.NormalizerInput

combined with other characters to make new characters. For example, the letter E combines with an acute accent to make E-acute (É). In some cases, the combined form (E-acute) should be treated as equivalent to the - uncombined form (E plus acute). For this reason, the standard defines four + uncombined form (E plus acute). Therefore, the standard defines four normalization forms that convert strings to a single equivalent form:

PreRound

-
public EDecimal PreRound(EContext ctx) Returns a number in which the value of this object is rounded to fit the  maximum precision allowed if it has more significant digits than the maximum  precision. The maximum precision allowed is given in an arithmetic context.  This method is designed for preparing operands to a custom arithmetic  operation in accordance with the "simplified" arithmetic given in Appendix A  of the General Decimal Arithmetic Specification.
+
public EDecimal PreRound(EContext ctx) Returns a number in which the value of this object is rounded to fit the  maximum precision allowed if it has more significant digits than the maximum  precision. The maximum precision allowed is given in an arithmetic context.  This method is designed for preparing operands to a custom arithmetic  operation per the "simplified" arithmetic given in Appendix A  of the General Decimal Arithmetic Specification.
 

Parameters:

diff --git a/Numbers/api/com.upokecenter.numbers.EFloat.html b/Numbers/api/com.upokecenter.numbers.EFloat.html index 5646b75cb..73601a716 100755 --- a/Numbers/api/com.upokecenter.numbers.EFloat.html +++ b/Numbers/api/com.upokecenter.numbers.EFloat.html @@ -4916,7 +4916,7 @@

RoundToPrecision

PreRound

-
public EFloat PreRound(EContext ctx) Returns a number in which the value of this object is rounded to fit the  maximum precision allowed if it has more significant digits than the maximum  precision. The maximum precision allowed is given in an arithmetic context.  This method is designed for preparing operands to a custom arithmetic  operation in accordance with the "simplified" arithmetic given in Appendix A  of the General Decimal Arithmetic Specification.
+
public EFloat PreRound(EContext ctx) Returns a number in which the value of this object is rounded to fit the  maximum precision allowed if it has more significant digits than the maximum  precision. The maximum precision allowed is given in an arithmetic context.  This method is designed for preparing operands to a custom arithmetic  operation per the "simplified" arithmetic given in Appendix A  of the General Decimal Arithmetic Specification.
 

Parameters:

diff --git a/Numbers/api/com.upokecenter.numbers.ERational.html b/Numbers/api/com.upokecenter.numbers.ERational.html index 1d6a9855e..e51ef1944 100755 --- a/Numbers/api/com.upokecenter.numbers.ERational.html +++ b/Numbers/api/com.upokecenter.numbers.ERational.html @@ -1236,7 +1236,7 @@

FromString

sequence of char s can also be "-INF", "-Infinity", "Infinity", "INF", quiet NaN ("NaN" /"-NaN") followed by any number of digits, or signaling NaN ("sNaN" /"-sNaN") followed by any number of digits, all in any - combination of upper and lowercase.

+ combination of uppercase and lowercase.

All characters mentioned earlier are the corresponding characters in the Basic Latin range. In particular, the digits must be the basic digits 0 to 9 (U+0030 to U+0039). The sequence @@ -1323,8 +1323,7 @@

FromString

The sequence of bytes can also be "-INF", "-Infinity", "Infinity", "INF", quiet NaN ("NaN" /"-NaN") followed by any number of digits, or signaling NaN ("sNaN" - /"-sNaN") followed by any number of digits, all in any combination of upper - and lower case.

+ /"-sNaN") followed by any number of digits, all in any combination of uppercase and lowercase.

All characters mentioned earlier are the corresponding characters in the Basic Latin range. In particular, the digits must be the basic digits 0 to 9 (U+0030 to U+0039). The sequence of bytes is diff --git a/betadist.py b/betadist.py index fccee5121..4e4d26c73 100755 --- a/betadist.py +++ b/betadist.py @@ -14,7 +14,6 @@ import math from fractions import Fraction - def betabin(k, psi, rho, cpsi, m=5): ret = math.comb(m - 1, k - 1) for i in range(k - 1): @@ -23,7 +22,6 @@ def betabin(k, psi, rho, cpsi, m=5): ret *= (m - psi) * rho / (m - 1) + j * (cpsi - rho) return ret - def genscore_mean_var(mean, vari, m=5): # Generalized score distribution, parameterized by mean and variance psi = mean @@ -36,7 +34,6 @@ def genscore_mean_var(mean, vari, m=5): else: return genscore(psi, 1 - (vari - vmin) / (vmax - vmin), m=m) - def genscore(psi, rho, m=5): # Generalized score distribution (GSD). # Return value is an integer in [1, m]. @@ -84,7 +81,6 @@ def genscore(psi, rho, m=5): else: return 1 + sum(randBernoulli((psi - 1) / (m - 1)) for i in range(m - 1)) - def randBernoulli(f): if f == 1: return 1 @@ -96,7 +92,6 @@ def randBernoulli(f): return 1 if realIsLess(RandUniform(), f) else 0 return 1 if realIsLess(RandUniform(), RealFraction(f)) else 0 - def tulap(m, b, q): # Tulap(m, b, q). ("truncated uniform Laplace"), m real, b in (0, 1), q in [0, 1). # Discrete Laplace(b) is Tulap(0,b,0) rounded to nearest integer. @@ -133,7 +128,6 @@ def tulap(m, b, q): if realIsLess(q / 2, fn) and realIsLess(fn, 1 - q / 2): return nreal + m - def gen_to_transition(s): size = len(s) m = [[0 for i in range(size + 1)] for j in range(size)] @@ -146,7 +140,6 @@ def gen_to_transition(s): m[i][j] = -s[i][j] / s[i][i] return m - class PhaseType: # Samples from a continuous phase-type distribution. # N is the number of states (other than the @@ -168,7 +161,6 @@ def sample(self): state = random.choices(range(self.n + 1), weights=self.trans[state])[0] return ret - def exchangeable_bernoulli(p, d, lamda=None): # p=expected value (in [0, 1]); d=dimension; lamda=weights for # each ray density @@ -209,10 +201,8 @@ def exchangeable_bernoulli(p, d, lamda=None): random.shuffle(ret) return ret - ################### - def psrn_complement(x): # NOTE: Assumes digits is 2 for i in range(len(x[2])): @@ -220,11 +210,9 @@ def psrn_complement(x): x[2][i] = 1 - x[2][i] return x - def psrn_new_01(): return [1, 0, []] - def psrn_fill(rg, psrn, precision=53, digits=2): af = 0 afrac = psrn[2] @@ -264,7 +252,6 @@ def psrn_fill(rg, psrn, precision=53, digits=2): asign * ((af // digits) + (aint * digits**precision)) / (digits**precision) ) - def psrn_in_range(rg, bmin, bmax, digits=2): if bmin >= bmax: raise ValueError @@ -304,7 +291,6 @@ def psrn_in_range(rg, bmin, bmax, digits=2): a[0] = 1 return a - def psrn_in_range_positive(rg, bmin, bmax, digits=2): if bmin >= bmax or bmin < 0 or bmax <= 0: raise ValueError @@ -404,11 +390,9 @@ def psrn_in_range_positive(rg, bmin, bmax, digits=2): continue return a - def psrn_sample(rg, psrn, digits=2): return psrn_less(rg, psrn_new_01(), [1, 0, psrn[2]], digits) - def psrn_less(rg, psrn1, psrn2, digits=2): if psrn1[0] == None or psrn1[1] == None or psrn2[0] == None or psrn2[1] == None: raise ValueError @@ -448,7 +432,6 @@ def psrn_less(rg, psrn1, psrn2, digits=2): return 0 index += 1 - def psrn_less_than_fraction(rg, psrn, rat, digits=2): if (psrn[0] != -1 and psrn[0] != 1) or psrn[1] == None: raise ValueError @@ -512,7 +495,6 @@ def psrn_less_than_fraction(rg, psrn, rat, digits=2): pt *= digits index += 1 - def psrn_reciprocal(rg, psrn1, digits=2): """Generates the reciprocal of a partially-sampled random number. psrn1: List containing the sign, integer part, and fractional part @@ -581,7 +563,6 @@ def psrn_reciprocal(rg, psrn1, digits=2): ddc *= digits dc *= digits - def proddist(x, a, b, c, d): if a * d < b * c: aa = a @@ -604,7 +585,6 @@ def proddist(x, a, b, c, d): r = max(0, min(1, math.log(b * d / x))) / math.log(b / a) return max(0, min(1, r)) - def proddist2(x, a, b, c, d): if a * d < b * c: aa = a @@ -627,7 +607,6 @@ def proddist2(x, a, b, c, d): r = [Fraction(b * d, x), Fraction(b, a)] return r - def psrn_multiply(rg, psrn1, psrn2, digits=2): """Multiplies two uniform partially-sampled random numbers. psrn1: List containing the sign, integer part, and fractional part @@ -638,7 +617,6 @@ def psrn_multiply(rg, psrn1, psrn2, digits=2): digits: Digit base of PSRNs' digits. Default is 2, or binary.""" return psrn_multiply_b(rg, psrn1, psrn2, digits=digits) - def _dlc(rg, psrn, c, digits=2): i = rg.rndint(c - 1) if i < psrn[1]: @@ -647,7 +625,6 @@ def _dlc(rg, psrn, c, digits=2): return psrn_sample(rg, psrn, digits=digits) return 0 - def _log_1n(rg, c): # ln(1+1/c) u = None @@ -660,7 +637,6 @@ def _log_1n(rg, c): if rg.rndint(c - 1) == 0: return 0 - def _log_yxyz(rg, psrn, st, digits=2): # ln((st+psrn)/st) if psrn[0] < 0 or psrn[1] > st: @@ -675,7 +651,6 @@ def _log_yxyz(rg, psrn, st, digits=2): if _dlc(rg, psrn, st, digits=digits) == 1: return 0 - def _dlc2(rg, psrn, n, d, digits=2): while True: if rg.rndint(d) != d: @@ -688,7 +663,6 @@ def _dlc2(rg, psrn, n, d, digits=2): if psrn_sample(rg, psrn, digits=digits) == 1: return 0 - def _log_xyzy(rg, psrn, large, midmax, digits=2): # ln(large/(midmax+psrn)) if psrn[0] < 0: @@ -705,7 +679,6 @@ def _log_xyzy(rg, psrn, large, midmax, digits=2): if _dlc2(rg, psrn, n, d, digits=digits) == 1: return 0 - def _log_yxyz_test(ps, st): import scipy.integrate as spi @@ -713,7 +686,6 @@ def _log_yxyz_test(ps, st): print(spi.quad(lambda x: math.log((st + ps + x) / st), 0, 1)[0]) print(sum(_log_yxyz(rg, [1, ps, []], st) for i in range(100000)) / 100000) - def _log_xyzy_test(ps, large, midmax): import scipy.integrate as spi @@ -723,7 +695,6 @@ def _log_xyzy_test(ps, large, midmax): sum(_log_xyzy(rg, [1, ps, []], large, midmax) for i in range(100000)) / 100000 ) - def _dlc2_test(n, d): import scipy.integrate as spi @@ -731,7 +702,6 @@ def _dlc2_test(n, d): print(spi.quad(lambda x: (n + (1 - x)) / (d + x), 0, 1)[0]) print(sum(_dlc2(rg, [1, 0, []], n, d) for i in range(100000)) / 100000) - if False: _dlc2_test(0, 19) _dlc2_test(12, 19) @@ -747,7 +717,6 @@ def _dlc2_test(n, d): _log_xyzy_test(0, 16, 14) exit() - def psrn_multiply_b(rg, psrn1, psrn2, digits=2, testing=False): if psrn1[0] == None or psrn1[1] == None or psrn2[0] == None or psrn2[1] == None: raise ValueError @@ -839,7 +808,6 @@ def psrn_multiply_b(rg, psrn1, psrn2, digits=2, testing=False): # if iters>100:print(iters) return cpsrn - def psrn_multiply_by_fraction(rg, psrn1, fraction, digits=2): """Multiplies a partially-sampled random number by a fraction. psrn1: List containing the sign, integer part, and fractional part @@ -886,7 +854,6 @@ def psrn_multiply_by_fraction(rg, psrn1, fraction, digits=2): dcount += 1 ddc *= digits - def psrn_add(rg, psrn1, psrn2, digits=2): """Adds two uniform partially-sampled random numbers. psrn1: List containing the sign, integer part, and fractional part @@ -989,7 +956,6 @@ def psrn_add(rg, psrn1, psrn2, digits=2): b *= digits newdigits += 1 - def psrn_add_fraction(rg, psrn, fraction, digits=2): if psrn[0] == None or psrn[1] == None: raise ValueError @@ -1098,7 +1064,6 @@ def psrn_add_fraction(rg, psrn, fraction, digits=2): rv = rv * digits + rg.rndint(digits - 1) rvs = rv + rvstart - def psrnexpo(rg): count = 0 while True: @@ -1116,10 +1081,8 @@ def psrnexpo(rg): return [1, count, y1[2]] count += 1 - ################### - def geobagcompare(bag, f): """Returns 1 with probability f(U), where U is the value that the given geometric bag turns out to hold, or 0 otherwise. @@ -1189,7 +1152,6 @@ def geobagcompare(bag, f): i += 1 iprec = 1 << i - def _bern_power(bern, bag, num, den, bagfactory): if len(bag) >= 4 and bag[0] == 0 and bag[1] == 0 and bag[2] == 0 and bag[3] == 0: # If the geometric bag is known to hold a very small number, use @@ -1203,14 +1165,12 @@ def _bern_power(bern, bag, num, den, bagfactory): else: return bern.power(bagfactory, num, den) - def _power_of_uniform_greaterthan1(bern, power, complement=False, precision=53): return psrn_fill( _power_of_uniform_greaterthan1_geobag(bern, power, complement), precision=precision, ) - def _power_of_uniform_greaterthan1_geobag(bern, power, complement=False): if power < 1: raise ValueError("Not supported") @@ -1250,7 +1210,6 @@ def _power_of_uniform_greaterthan1_geobag(bern, power, complement=False): # Flip all bits if complement is true return psrn_complement(bag) if complement else bag - def powerOfUniform(b, px, py, precision=53): """Generates a power of a uniform random number. - px, py - Numerator and denominator of desired exponent for the uniform @@ -1262,11 +1221,9 @@ def powerOfUniform(b, px, py, precision=53): # is in (0, 1]. return betadist(b, py, px, 1, 1, precision) - def betadist(b, ax=1, ay=1, bx=1, by=1, precision=53): return psrn_fill(betadist_geobag(b, ax, ay, bx, by), precision=precision) - def betadist_geobag(b, ax=1, ay=1, bx=1, by=1): """Generates a beta-distributed random number with arbitrary (user-defined) precision. Currently, this sampler only works if (ax/ay) and @@ -1329,10 +1286,8 @@ def betadist_geobag(b, ax=1, ay=1, bx=1, by=1): # Accepted return bag - ##################### - class _RGConv: # Wrapper that takes an object that implements # rndint(maxinc). Implements randint(mininc,maxinc). @@ -1342,7 +1297,6 @@ def __init__(self, rg): def randint(self, a, b): return a + self.rg.rndint(b - a) - ##################### def truncated_gamma(rg, bern, ax, ay, precision=53): # Văduva's gamma generator truncated to [0, 1], @@ -1360,7 +1314,6 @@ def truncated_gamma(rg, bern, ax, ay, precision=53): w = u k += 1 - def forsythe_prob2(rg, x): # Returns 1 with probability x*exp(1-x), where x is in [0, 1]. # Implemented with the help of Theorem IV.2.1(iii) given in @@ -1385,7 +1338,6 @@ def forsythe_prob2(rg, x): if k % 2 == 1: return 1 if psrn_less_than_fraction(rg, ret, x) == 1 else 0 - def forsythe_prob3(rg, x): # Returns 1 with probability erf(x)/erf(1), where x is in [0, 1]. # Implemented with the help of Theorem IV.2.1(iii) given in @@ -1407,7 +1359,6 @@ def forsythe_prob3(rg, x): if k % 2 == 1: return 1 if psrn_less_than_fraction(rg, ret, x) == 1 else 0 - def forsythe_prob(rg, m, n): # Returns 1 with probability gamma(m,n)/gamma(m,1), # where gamma(.) is the lower incomplete gamma function. @@ -1428,7 +1379,6 @@ def forsythe_prob(rg, m, n): if k % 2 == 1: return 1 if psrn_less_than_fraction(rg, ret, n) == 1 else 0 - def exp_minus_x2y(rg, f, y, pwr=2): """B(x) -> B(exp(-x*x*y))""" if y > 1: @@ -1463,7 +1413,6 @@ def exp_minus_x2y(rg, f, y, pwr=2): fac *= n y *= uy - def exp_minus_xy(rg, f, y): """B(x) -> B(exp(-x*y))""" if y > 1: @@ -1497,7 +1446,6 @@ def exp_minus_xy(rg, f, y): fac *= n y *= uy - def sampleIntPlusBag(rg, psrn, k): """Return 1 with probability (x+k)/2^bitlength(k). Ignores PSRN's integer part and sign.""" @@ -1519,7 +1467,6 @@ def sampleIntPlusBag(rg, psrn, k): psrn[2][r] = rg.randbit() return psrn[2][r] - def rayleighpsrn(rg, s=1): k = 0 y = Fraction(2 * s * s) @@ -1552,7 +1499,6 @@ def rayleighpsrn(rg, s=1): bag[1] = k return bag - def size_biased_poisson_ailamujia(rg, eta=1): """Hassan, A., Dar, S.A., et al., "On size biased Poisson Ailamujia distribution and its applications", Pak. J. Statistics 37(1), 19-38, 2021.""" @@ -1568,7 +1514,6 @@ def size_biased_poisson_ailamujia(rg, eta=1): cumu -= pr z += 1 - def genshape(rg, inshape): """Generates a random point inside a 2-dimensional shape, in the form of a uniform PSRN. inshape is a function that takes three parameters (x, y, s) and @@ -1603,7 +1548,6 @@ def genshape(rg, inshape): s *= base d += 1 - class _PavingNode: def __init__(self, coords): """Represents a tree node of a mapped regular paving.""" @@ -1640,7 +1584,6 @@ def bisect(self): self.left, self.right = self.getBisection() return self.left, self.right - class ShapeSampler2: YES = 1 NO = -1 @@ -1730,7 +1673,6 @@ def sample(self, rg): box = left if child == 0 else right box.mark = self.inshape(box.coords) - class ShapeSampler: def __init__(self, inshape, dx=1, dy=1): """Builds a sampler for random numbers (in the form of PSRNs) on or inside a 2-dimensional shape. @@ -1803,7 +1745,6 @@ def sample(self, rg): psrny[0] = -1 if rg.randbit() == 0 else 1 return [psrnx, psrny] - def _peres(bits, output): u = [] v = [] @@ -1829,7 +1770,6 @@ def _peres(bits, output): _peres(u, output) _peres(v, output) - class _BitFetchingRandomGen: def __init__(self, *args): self.rg = randomgen.RandomGen(*args) @@ -1907,7 +1847,6 @@ def feed_fetchedbits(self): x >>= 1 self.fetchedbits = 0 - def _test_rand_extraction(rg, func, digits=2, nofill=False): func(rg) return @@ -1945,7 +1884,6 @@ def _test_rand_extraction(rg, func, digits=2, nofill=False): except: pass - def addto1(rg): # Simulates the number of uniforms # in [0,1] needed to produce a sum greater than 1. @@ -1958,10 +1896,8 @@ def addto1(rg): break return i - ######################## - def recordcount(n): if n == 1: return 1 @@ -1971,7 +1907,6 @@ def recordcount(n): ret += 1 return ret - def lah(n, k): # Lah distribution. # Kabluchko, Zakhar, and Alexander Marynych. @@ -2000,7 +1935,6 @@ def lah(n, k): run += 1 return ret - # Random real numbers. Uses ideas from H.J. Boehm's # concept of "constructive reals". # Essentially, each operation on real numbers has @@ -2057,7 +1991,6 @@ def lah(n, k): # shared with other tasks). # - class Real: def ev(self, n): raise NotImplementedError @@ -2149,7 +2082,6 @@ def evstable(a, prec): def __repr__(a): return "Real" - class RealPi(Real): def __init__(self, fraction=1, consistent=False): self.fraction = Fraction(fraction) @@ -2207,11 +2139,9 @@ def ev(self, n): lower = upper k += 1 - REALPI = RealPi(consistent=True) REALHALFPI = RealPi(Fraction(1, 2), consistent=True) - class RealTan(Real): def __init__(self, a): self.a = a if isinstance(a, Real) else RealFraction(a) @@ -2225,7 +2155,6 @@ def __repr__(self): def ev(self, n): return self.a.ev(n) - class RealArcTan2(Real): def __init__(self, y, x): self.r = None @@ -2258,7 +2187,6 @@ def __repr__(self): def ev(self, n): return self.r.ev(n) - _BERNNUMBERS = [ Fraction(1), Fraction(1, 2), @@ -2275,7 +2203,6 @@ def ev(self, n): ] _extrabernnumbers = {} - def bernoullinum(n): # Calculates Bernoulli numbers global _extrabernnumbers @@ -2296,10 +2223,8 @@ def bernoullinum(n): _extrabernnumbers[n] = ret return ret - _STIRLING1 = {} - def stirling1(n, k): # Calculates Stirling numbers of the first kind if n == k: @@ -2312,7 +2237,6 @@ def stirling1(n, k): _STIRLING1[(n, k)] = ret return ret - def loggammahelper(n, precision): # Implements the first listed convergent version of Stirling's formula # given in section 3 of R. Schumacher, "Rapidly Convergent Summation Formulas @@ -2373,7 +2297,6 @@ def loggammahelper(n, precision): errintv.addnumden(result.supn, result.supd) return (result.infn, result.infd, errintv.supn, errintv.supd) - class RealLogGammaInt(Real): _logPi = None @@ -2406,7 +2329,6 @@ def ev(self, n): return 0 return self.r.ev(n) - class _RealLogGammaIntHelper(Real): def __init__(self, a): self.a = a @@ -2447,7 +2369,6 @@ def ev(self, n): return cinf nv += 6 - def logbinco(n, k): # Log binomial coefficient. if k + 1 == (n - k) + 1: @@ -2460,19 +2381,16 @@ def logbinco(n, k): ) return r - def logbinprob(n, k): # Log of binomial probability, that is, the log of the probability # that exactly k zeros occur among n unbiased random bits. divisor = RealLn(2) * n # ln(2)*n = ln(2**n) return logbinco(n, k) - divisor - def logpoisson(lamda, n): # Log of the probability that a Poisson(lamda) random number is n. return RealLn(lamda) * n - lamda - RealLogGammaInt(n + 1) - class RealErf(Real): def __init__(self, a): self.a = a @@ -2484,7 +2402,6 @@ def __repr__(self): def ev(self, n): return self.r.ev(n) - class _RealErfHelper(Real): def __init__(self, a): self.r = None @@ -2582,7 +2499,6 @@ def ev(self, n): return cinf nv += 6 - class RealArcTan(Real): def __init__(self, a): self.r = None @@ -2700,7 +2616,6 @@ def ev(self, n): return cinf nv += 6 - class RealSin(Real): def __init__(self, a): self.a = a if isinstance(a, Real) else RealFraction(a) @@ -2720,7 +2635,6 @@ def ev(self, n): return 0 return self.c.ev(n) - def fracAreClose(a, b, n): # if a>b: raise ValueError an = a.numerator @@ -2729,7 +2643,6 @@ def fracAreClose(a, b, n): bd = b.denominator return fracAreCloseND(an, ad, bn, bd, n) - def fracEV(sn, sd, n): sn = sn << (n + 2) if sn >= 0 else sn * (1 << (n + 2)) ret = abs(sn) // sd @@ -2743,7 +2656,6 @@ def fracEV(sn, sd, n): ret2 = (ret // 4) + 1 if ret % 4 >= 2 else (ret // 4) return ret2 - def fracAreCloseND(an, ad, bn, bd, n): if ad < 0 or bd < 0: raise ValueError @@ -2766,7 +2678,6 @@ def fracAreCloseND(an, ad, bn, bd, n): return ra return None - class RealCos(Real): def __init__(self, a): # NOTE: Works best if 'a' is in (-2*pi, 2*pi). @@ -2870,7 +2781,6 @@ def ev(self, n): return cinf nv += 6 - class RealExp(Real): def __init__(self, a): self.a = a if isinstance(a, Real) else RealFraction(a) @@ -2976,7 +2886,6 @@ def ev(self, n): return cinf nv += 2 - def randMax(n=2): # rand()^(1/n) u = RandUniform() for i in range(1, n): @@ -2985,7 +2894,6 @@ def randMax(n=2): # rand()^(1/n) u = v return u - def randMin(n=2): # 1-rand()^(1/n) u = RandUniform() for i in range(1, n): @@ -2994,7 +2902,6 @@ def randMin(n=2): # 1-rand()^(1/n) u = v return u - class RealPow(Real): def __init__(self, a, b): self.powint = None @@ -3130,7 +3037,6 @@ def ev(self, n): nv += 4 return self.r.ev(n) - class RealDivide(Real): def __init__(self, a, b): self.a = a if isinstance(a, Real) else RealFraction(a) @@ -3208,7 +3114,6 @@ def ev(self, n): return cinf nv += 4 - class RandPSRN(Real): # Constructive real operation # that takes a positive uniform PSRN with base-2 @@ -3236,10 +3141,8 @@ def ev(self, n): bits = -bits return bits - _realbits = 0 - class RandUniform(Real): # Random uniform real number in the interval (0, 1). def __init__(self): @@ -3267,7 +3170,6 @@ def ev(self, n): # to strictly less than 1 ulp return ret - class RealFraction(Real): def __init__(self, a, b=None): if isinstance(a, Real): @@ -3401,7 +3303,6 @@ def ev(self, n): return self.num << n return fracEV(self.num, self.den, n) - class RealNegate(Real): def __init__(self, a): self.a = a if isinstance(a, Real) else RealFraction(a) @@ -3412,7 +3313,6 @@ def __repr__(self): def ev(self, n): return -self.a.ev(n) - class RealSubtract(Real): def __init__(self, a, b): self.a = a if isinstance(a, Real) else RealFraction(a) @@ -3425,7 +3325,6 @@ def ev(self, n): r = self.a.ev(n + 2) - self.b.ev(n + 2) return (r // 4) + 1 if r % 4 >= 2 else (r // 4) - class RealAdd(Real): def __init__(self, a, b): self.a = a if isinstance(a, Real) else RealFraction(a) @@ -3455,7 +3354,6 @@ def ev(self, n): self.ev_v = ret return ret - # floor(atanh(1/2^i)*2^29) ArcTanHTable = [ 0, # Infinity @@ -3490,7 +3388,6 @@ def ev(self, n): 1, ] - def _roundedshiftraw(a, shift): aa = abs(a) ret = aa >> shift @@ -3502,7 +3399,6 @@ def _roundedshiftraw(a, shift): ret = -ret return ret - def logsmall(av, n): if n > 16 or av < 2: return None @@ -3518,14 +3414,12 @@ def logsmall(av, n): supcr = supcr + 2 if avplus <= 271 else supcr + 1 return (infcr, 65536, supcr, 65536) - CRUDELOG_BITS = 16 CRUDELOG_LOGMIN = (1 << CRUDELOG_BITS) * 15 // 100 CRUDELOG_LOG2BITS = 45426 CRUDELOG_ARCTANFRAC = 29 CRUDELOG_ARCTANBITDIFF = CRUDELOG_ARCTANFRAC - CRUDELOG_BITS - def crudelog(av): # The CORDIC way to # calculate an approximation of the natural logarithm of av/2^_bits. @@ -3559,7 +3453,6 @@ def crudelog(av): rz += ArcTanHTable[i] return _roundedshiftraw(rz, CRUDELOG_ARCTANBITDIFF - 1) - CRUDELOG = [0 if i == 0 else crudelog(i) for i in range(0, 65536 + 1)] # 10-degree Chebyshev approximation of ln(x) on [1, 1.0625]. @@ -3660,7 +3553,6 @@ def crudelog(av): ), ] - class FPInterval: # Fixed-point interval def __init__(self, n, d, prec): @@ -3782,7 +3674,6 @@ def mulnumden(self, n, d): self._truncate() # print(["mul1",self.infn/self.infd,self.supn/self.supd]) - class RealLn(Real): def __init__(self, a): self.a = a if isinstance(a, Real) else RealFraction(a) @@ -3987,7 +3878,6 @@ def ev(self, n): return cinf nv += 4 - def realFloor(a): if isinstance(a, RealFraction): # NOTE: Python's "//" operator does floor division @@ -4004,14 +3894,12 @@ def realFloor(a): return av // (1 << n) n += 2 - def realCeiling(a): if isinstance(a, RealFraction): # NOTE: Python's "//" operator does floor division return -((-a.num) // a.den) return -realFloor(-a) - def realIsLessOrEqual(a, b): a = a if isinstance(a, Real) else RealFraction(a) b = b if isinstance(b, Real) else RealFraction(b) @@ -4019,7 +3907,6 @@ def realIsLessOrEqual(a, b): return Fraction(a.num, a.den) <= Fraction(b.num, b.den) return realIsLess(a, b) - def realIsLess(a, b): if isinstance(a, int): a = Fraction(a) @@ -4049,11 +3936,9 @@ def realIsLess(a, b): return True n += 3 - def realIsGreater(a, b): return realIsLess(b, a) - def realIsNegative(a): try: return a.isNegative() @@ -4069,7 +3954,6 @@ def realIsNegative(a): return True n += 2 - class RealSqrt(Real): def __init__(self, a): self.a = RealPow(a, Fraction(1, 2)) @@ -4080,7 +3964,6 @@ def __repr__(self): def ev(self, n): return self.a.ev(n) - class RealMultiply(Real): def __init__(self, a, b): self.a = a if isinstance(a, Real) else RealFraction(a) @@ -4191,10 +4074,8 @@ def ev(self, n): return cinf nv += 6 - REAL_858_1000 = RealFraction(Fraction(858, 1000)) - def realNormalROU(mu=0, sigma=1): # Generates a Gaussian random variate using # the ratio of uniforms method. @@ -4209,7 +4090,6 @@ def realNormalROU(mu=0, sigma=1): return -b / a if mu == 0 else (-b / a) + mu return b / a if mu == 0 else (b / a) + mu - def fpNormalROU(mu=0, sigma=1): # Generates a Gaussian random variate using # the ratio of uniforms method. @@ -4226,7 +4106,6 @@ def fpNormalROU(mu=0, sigma=1): return -b / a if mu == 0 else (-b / a) + mu return b / a if mu == 0 else (b / a) + mu - def logconcave(f, c): # Devroye 1986, chapter 7 # Samples a random variate from an absolutely # continuous log-concave distribution. @@ -4245,7 +4124,6 @@ def logconcave(f, c): # Devroye 1986, chapter 7 if realIsLess(z, f(x)): return x - def gammaDist2(): # Gamma distribution with parameter 2, # using Devroye's algorithm for log concave densities @@ -4263,7 +4141,6 @@ def gammaDist2(): if realIsLess(zero, x) and z <= RealLn(x) - x + 1: return x - def c2a(r=None): # Generates a uniform random point on the unit circle # with radius r (default is 1). @@ -4284,7 +4161,6 @@ def c2a(r=None): else: [(y1 - y2) / s, 2 * x1 * x2 / s] - def muth(mu): # Ratio-of-uniforms sampler for the Muth distribution # Irshad, M. R., Maya, R., & Arun, S. P. (2021). Muth distribution @@ -4309,7 +4185,6 @@ def muth(mu): if realIsLess(RealLn(v) * 2, logpdf): return x - def realGamma(ml): # Generates a gamma random variate # using the Marsaglia--Tsang (2000) algorithm. @@ -4339,7 +4214,6 @@ def realGamma(ml): ret = ret * RealPow(RandUniform(), Fraction(1, ml)) return ret - class RandUniformIntFrac(Real): def __init__(self, i, f): if i < 0: @@ -4353,7 +4227,6 @@ def __repr__(self): def ev(self, n): return (self.i << n) + self.f.ev(n) - class RandUniformNegIntFrac(Real): def __init__(self, i, f): if i < 0: @@ -4367,7 +4240,6 @@ def __repr__(self): def ev(self, n): return -((self.i << n) + self.f.ev(n)) - def randLnUniform(): count = 0 while True: @@ -4387,7 +4259,6 @@ def randLnUniform(): return RandUniformNegIntFrac(count, y1) count += 1 - def randUniformLessThan(val): # NOTE: Assumes 'val' is greater than 0 tev = val.ev(0) + 1 @@ -4396,7 +4267,6 @@ def randUniformLessThan(val): if realIsLess(x, val): return x - def randUniformPower(pwr): if isinstance(pwr, RealFraction) and pwr.den > 0: if pwr.num < 0: @@ -4405,7 +4275,6 @@ def randUniformPower(pwr): return randMax(pwr.den) ** abs(pwr.num) return RandUniform() ** pwr - def monoSecondMoment(secondMoment, pdf): # For distributions on the positive real line # with: @@ -4443,7 +4312,6 @@ def monoSecondMoment(secondMoment, pdf): # print(iters) return x - class SinFunction: def value(self, pt): return ( @@ -4452,7 +4320,6 @@ def value(self, pt): else RealSin(pt) ) - def bernsteinDiff(coeffs, diff): # Gets the Bernstein coefficients of the 'diff'-th derivative of # a polynomial given its Bernstein coefficients. @@ -4472,7 +4339,6 @@ def bernsteinDiff(coeffs, diff): n -= 1 return coeffs - class BernsteinPoly: def __init__(self, coeffs): # Polynomial in Bernstein form defined on the closed unit interval. @@ -4523,7 +4389,6 @@ def value(self, pt): def diff(self, pt, d=1): return self.deriv(d).value(pt) - def minDegree(maxValue, maxDeriv, epsilon, deriv=4): # Minimum degree for iteratedPoly2 (for deriv=4) or # iteratedPoly3 (for deriv=5 or 6) needed to achieve @@ -4553,7 +4418,6 @@ def minDegree(maxValue, maxDeriv, epsilon, deriv=4): ) raise ValueError("'deriv' value not supported") - def iteratedPoly2(func, n): # Bernstein coefficients for the # Micchelli-Felbecker iterated Bernstein polynomial of order 2 @@ -4564,7 +4428,6 @@ def iteratedPoly2(func, n): ret.append(func.value(rf) * 2 - bp.value(rf)) return ret - def iteratedPoly3(func, n): # Bernstein coefficients for the # Micchelli-Felbecker iterated Bernstein polynomial of order 3 @@ -4576,7 +4439,6 @@ def iteratedPoly3(func, n): ret.append(bprec.value(rf) + 3 * (func.value(rf) - bp.value(rf))) return ret - class PiecewiseBernstein: def __init__(self): self.pieces = [] @@ -4640,7 +4502,6 @@ def get_coeffs(self): raise ValueError("likely not a polynomial") return [v for v in self.pieces[0][0]] - def c4example(): # Example function: A C4 continuous piecewise polynomial return ( @@ -4670,7 +4531,6 @@ def c4example(): ) ) - def iteratedPolyExample(): # Example of a Bernoulli factory that simulates # a polynomial that is close to a C4 continuous function @@ -4688,7 +4548,6 @@ def iteratedPolyExample(): print(r / 10000) print(sinf.value(RealFraction(0.3)).disp()) - ###################### if __name__ == "__main__": diff --git a/binomial.py b/binomial.py index f1179dbc9..e03be0e3a 100755 --- a/binomial.py +++ b/binomial.py @@ -13,7 +13,6 @@ realIsLess, ) - class BinomialSampler: def __init__(self, rg=None): self.rg = randomgen.RandomGen() if rg == None else rg diff --git a/exporand.html b/exporand.html index cb5759027..742114d86 100755 --- a/exporand.html +++ b/exporand.html @@ -818,7 +818,7 @@

kthsmallest

-

Power-of-Uniform Sub-Algorithm

+

Power-of-Uniform Sub-Algorithm

The power-of-uniform subalgorithm is used for certain cases of the beta sampler below. It returns Upx/py, where U is a uniform random variate in the interval [0, 1] and px/py is greater than 1, but unlike the naïve algorithm it supports an arbitrary precision, uses only random bits, and avoids floating-point arithmetic. It also uses a complement flag to determine whether to return 1 minus the result.

@@ -3306,7 +3306,7 @@

License

Hill, T.P. and Schürger, K., 2005. Regularity of digits and significant digits of random variables. Stochastic processes and their applications, 115(10), pp.1723-1743. 

  • -

    A non-discrete distribution is a probability distribution taking on values that each can’t be mapped to a different integer. An example is a distribution taking on any real number between 0 and 1. 

    +

    A nondiscrete distribution is a probability distribution taking on values that each can’t be mapped to a different integer. An example is a distribution taking on any real number between 0 and 1. 

  • J.F. Williamson, “Random selection of points distributed on curved surfaces”, Physics in Medicine & Biology 32(10), 1987. 

    diff --git a/fixed.py b/fixed.py index 17c265f6b..3b11c339e 100755 --- a/fixed.py +++ b/fixed.py @@ -1,7 +1,6 @@ import decimal as dec import math - class Fixed: """ Fixed-point numbers, represented using integers that store multiples @@ -758,7 +757,6 @@ def __str__(self): def __repr__(self): return str(dec.Decimal(self.value) / dec.Decimal((1 << Fixed.BITS))) - if __name__ == "__main__": def asserteq(a, b): diff --git a/html3dutil/extras_contourlines.html b/html3dutil/extras_contourlines.html index 839219dba..5b3a7eaa0 100755 --- a/html3dutil/extras_contourlines.html +++ b/html3dutil/extras_contourlines.html @@ -37,7 +37,7 @@

    Methods

    Parameters

      -
    • func (Type: function)
      A function that takes two parameters–a u-coordinate and a V coordinate–and returns a number at that point.
    • +
    • func (Type: function)
      A function that takes two parameters–a u-coordinate and a v-coordinate–and returns a number at that point.
    • levels (Type: Array.<number>)
      An array of values at which to draw contour lines. For example, if levels is [20, 25], this function will draw contour lines along the values 20 and 25.
    • u1 (Type: number)
      Starting u-coordinate to sample.
    • u2 (Type: number)
      Ending u-coordinate to sample.
    • diff --git a/interval.py b/interval.py index 55993abb5..cc662a51d 100755 --- a/interval.py +++ b/interval.py @@ -18,7 +18,6 @@ # precision in bits. # - class FInterval: """An interval of two Fractions. x.sup holds the upper bound, and x.inf holds the lower bound.""" @@ -344,7 +343,6 @@ def cos(self, precision): def __repr__(self): return "[%s, %s]" % (float(self.inf), float(self.sup)) - def _polynomialProduct(a, b): # Finds the product of two polynomials. Each polynomial # is a list of the following form: @@ -356,7 +354,6 @@ def _polynomialProduct(a, b): ret[i + j] += a[i] * b[j] return ret - def _polynomialProductA(a, b0, b1): ret = [0 for i in range(len(a) + 1)] for i in range(len(a)): @@ -364,7 +361,6 @@ def _polynomialProductA(a, b0, b1): ret[i + 1] += a[i] * b1 # b1 is 1st-order coefficient return ret - def _polynomialIntegral(p, x=1): # Finds the integral of a polynomial at the point x. if x == 1: diff --git a/keaneobrien.py b/keaneobrien.py index 676796a9c..cf9f05ba4 100755 --- a/keaneobrien.py +++ b/keaneobrien.py @@ -1,7 +1,6 @@ from sympy import * import random - def roughmaxiopen(func, x): # Rough version of maxiopen that samples points # in the open interval (0, 1) to find an approximate @@ -10,7 +9,6 @@ def roughmaxiopen(func, x): *[ceiling(func.subs(x, S(i) / 100) * 100000) / 100000 for i in range(1, 100)] ) - def maxi(func, x): # Finds the maximum of a function whose domain is [0,1]. # Unfortunately, SymPy's 'maximum' does not always work @@ -22,7 +20,6 @@ def maxi(func, x): # func=target function; x=variable used in 'func' return maximum(func, x, Interval(0, 1)) - def maxiopen(func, x): # Finds the maximum of a function whose domain is (0,1). # Unfortunately, SymPy's 'maximum' does not always work @@ -34,7 +31,6 @@ def maxiopen(func, x): # func=target function; x=variable used in 'func' return maximum(func, x, Interval.open(0, 1)) - def bernpoly(a, x): # Create a polynomial in Bernstein form using the given # array of coefficients and the symbol x. The polynomial's @@ -49,7 +45,6 @@ def bernpoly(a, x): ret += a[i] * bino * pt**i * (1 - pt) ** (n - i) return ret - def kobevents(func, p, numpolys=12): # Generates 'func', a Bernoulli factory function, # as a convex combination of polynomials @@ -150,7 +145,6 @@ def kobevents(func, p, numpolys=12): # plot(*[r[2] for r in results], (p, 0, 1)) return fevents - def kob(fevents, p): # Simulates a convex combination of polynomials # returned by kobevents, which will generally be an @@ -166,7 +160,6 @@ def kob(fevents, p): s = sum(1 if random.random() < p else 0 for i in range(len(fevents[geo]) - 1)) return fevents[geo][s] - p = symbols("p") # Variable used in the example function # Example function ofunc = sin(2 * pi * p) / 4 + S.Half + S(1) / 100 diff --git a/logconcave.py b/logconcave.py index 686000545..4542022e1 100755 --- a/logconcave.py +++ b/logconcave.py @@ -7,7 +7,6 @@ import math import random - class TConcaveDiscreteSampler: """ Generates a random variate that follows a @@ -100,7 +99,6 @@ def codegen(self, name, pdfcall=None): ret += " if u*u <= %s(ret): return ret\n\n" % (pdfcall) return ret - class TConcaveSampler: """ Generates a random variate that follows a @@ -232,7 +230,6 @@ def _simplesampleOne(self): if u * u <= self.f(ret): return ret - class LogConcaveSampler: """ Generates a random variate that follows a distribution @@ -357,7 +354,6 @@ def sampleOne(self): elif logchi <= pr: return ret - class UnimodalSampler: """ Generates a random variate that follows a unimodal distribution @@ -433,7 +429,6 @@ def sampleOne(self): self.mode, self.mode + 1, self.modecdf, self.modecdfright ) - class LogConcaveSampler2: """ Generates a random variate that follows a distribution @@ -490,7 +485,6 @@ def sampleOne(self): if ret != None: return ret - class LogConcaveSamplerMonotone: """ Generates a random variate that follows a distribution @@ -668,7 +662,6 @@ def sampleOne(self): if ret != None: return ret - class NormalDist: def __init__(self, mu=0, sigma=1): psi = lambda x: -(x * x) / (2 * sigma * sigma) @@ -678,7 +671,6 @@ def __init__(self, mu=0, sigma=1): def sample(self, n): return [x + self.mu for x in self.sampler.sample(n)] - class GenInvGaussianAlphaBeta: def __init__(self, lamda, alpha, beta): # Three-parameter version (Hörmann, W., Leydold, J., @@ -697,7 +689,6 @@ def fromLambdaPsiChi(lamda, psi, chi): lamda, math.sqrt(psi / chi), math.sqrt(psi * chi) ) - class GammaDist: # Devroye 2014 def __init__(self, a): @@ -722,7 +713,6 @@ def sample(self, n): for x in self.sampler.sample(n) ] - class GenInvGaussian: def __init__(self, lamda, omega): # omega>0, lamda>=0. Two-parameter version described in Devroye 2014 @@ -741,7 +731,6 @@ def _trans(self, v): def sample(self, n): return [self._trans(x) for x in self.sampler.sample(n)] - if __name__ == "__main__": def bucket(v, ls, buckets): diff --git a/lorentz.py b/lorentz.py index 6786739fb..3701d0f6b 100755 --- a/lorentz.py +++ b/lorentz.py @@ -7,7 +7,6 @@ BINCOMBS = {} SPLITBINCOMBS = {} - def ccomb(n, k): if k < 10: return math.comb(n, k) @@ -17,7 +16,6 @@ def ccomb(n, k): BINCOMBS[(n, k)] = v return v - def binco2(v, n, k): if k < 100 or k > n - 100: return v * ccomb(n, k) @@ -49,7 +47,6 @@ def binco2(v, n, k): v0 = v * s[1] return v0 - def polyshift(nrcoeffs, theta, d, alpha=2): # Upward and downward shift of polynomial according to step 5 # in Holtz et al. 2011, for even integer r>=2 or r=1 (r times @@ -79,7 +76,6 @@ def polyshift(nrcoeffs, theta, d, alpha=2): lower = [nrcoeffs[i] - phi[i] * d for i in range(len(phi))] return upper, lower - def example1(): # Example function: A concave piecewise # polynomial with three continuous derivatives @@ -107,11 +103,9 @@ def example1(): ) return pwp2 - REALONE = RealFraction(1) REALZERO = RealFraction(0) - def elevatemulti(coeffs, r): # Elevate polynomial in Bernstein form by r degrees if r < 0: raise ValueError @@ -127,7 +121,6 @@ def elevatemulti(coeffs, r): # Elevate polynomial in Bernstein form by r degree for k in range(n + r + 1) ] - def degelev(coeffs, r): # Elevate polynomial in Bernstein form by r degrees if r < 0: raise ValueError @@ -152,7 +145,6 @@ def degelev(coeffs, r): # Elevate polynomial in Bernstein form by r degrees ret[k] += coeffs[k - j] * binrj return [ret[i] / math.comb(n + r, i) for i in range(n + r + 1)] - class PolyDiff: def __init__(self, a, b): self.a = a @@ -164,7 +156,6 @@ def value(self, x): def diff(self, x, d=1): return self.a.diff(x, d=d) - self.b.diff(x, d=d) - def lorentz2iter(func, coeffs, nPlusR): existingNPlusR = len(coeffs) - 1 if existingNPlusR <= 0 or nPlusR < existingNPlusR: @@ -175,7 +166,6 @@ def lorentz2iter(func, coeffs, nPlusR): raise ValueError return [realSimplify(v + w) for v, w in zip(coeffs2, l4)] - def lorentz4iter(func, coeffs, nPlusR): existingNPlusR = len(coeffs) - 1 if existingNPlusR <= 0 or nPlusR < existingNPlusR: @@ -186,7 +176,6 @@ def lorentz4iter(func, coeffs, nPlusR): raise ValueError return [realSimplify(v + w) for v, w in zip(coeffs2, l4)] - def lorentz4poly(pwpoly, n): # Polynomial for Lorentz operator with r=4, # of degree n+r = n+4, given four times differentiable piecewise polynomial @@ -241,12 +230,10 @@ def lorentz4poly(pwpoly, n): # it into a Bernstein coefficient return [(vals[i]) / math.comb(n + r, i) for i in range(n + r + 1)] - # NOTE: Lorentz operator with r=0 or r=1 of degree n+r # is the same as the degree-n Bernstein polynomial, elevated # r degrees to degree n+r. - def lorentz2polyB(func, n, r=2): # Polynomial for Lorentz operator with r=0, 1, 2, # of degree n+r, given r times differentiable function @@ -269,7 +256,6 @@ def lorentz2polyB(func, n, r=2): vals[k] -= fv return vals - def lorentz2poly(pwpoly, n): # Polynomial for Lorentz operator with r=2, # of degree n+r = n+2, given twice differentiable piecewise polynomial @@ -307,7 +293,6 @@ def lorentz2poly(pwpoly, n): # it into a Bernstein coefficient return [(vals[i]) / math.comb(n + r, i) for i in range(n + r + 1)] - def realSimplify(r): if isinstance(r, RealMultiply): if isinstance(r.a, RealFraction) and isinstance(r.b, RealFraction): @@ -335,7 +320,6 @@ def realSimplify(r): return s1 - s2 return r - class C4Function: # func must map [0, 1] to (0, 1) and have at least # four continuous derivatives by default. @@ -520,13 +504,11 @@ def simulate(self, coin): """- coin(): Function that returns 1 or 0 with a fixed probability.""" return simulate(coin, self.fbelow, self.fabove, self.fbound, self.nextdegree) - def _fb2(fbelow, a, b, v=1): if b > a: raise ValueError return Fraction(realFloor(fbelow(a, b) * v), v) - def _fa2(fabove, a, b, v=1): if b > a: raise ValueError @@ -534,7 +516,6 @@ def _fa2(fabove, a, b, v=1): imv = int(mv) return Fraction(imv, v) if mv == imv else Fraction(imv + 1, v) - def simulate(coin, fbelow, fabove, fbound, nextdegree=None): """A Bernoulli factory for a continuous function f(x) that maps [0, 1] to [0, 1] (and where f(x) is polynomially bounded). Returns either 1 @@ -639,7 +620,6 @@ def simulate(coin, fbelow, fabove, fbound, nextdegree=None): lastdegree = degree degree = nextdegree(degree) if nextdegree != None else degree * 2 - def cc(): # ce = c4example() # f = C4Function(ce, 5, lorentz=False,contderivs=4) @@ -649,25 +629,21 @@ def cc(): print(ce.value(0.9)) print(sum(f.simulate(coin) for i in range(50000)) / 50000) - from sympy import Min, Max, ceiling, S from sympy import Matrix, binomial, chebyshevt, pi, Piecewise, Eq, floor, ceiling BOUNDMULT = 1000000000000000 - def upperbound(x, boundmult=BOUNDMULT): # Calculates a limited-precision upper bound of x. boundmult = S(boundmult) return S(int(ceiling(x * boundmult))) / boundmult - def lowerbound(x, boundmult=BOUNDMULT): # Calculates a limited-precision lower bound of x. boundmult = S(boundmult) return S(int(floor(x * boundmult))) / boundmult - class Intv: def __init__(self, v, w=None, boundmult=BOUNDMULT): if w != None: @@ -738,7 +714,6 @@ def __truediv__(a, b): a4 = a.sup * S(1) / b.sup return Intv(min(a1, a2, a3, a4), max(a1, a2, a3, a4)) - def tachevcoeffs(func, x, n): # Sample from func at (n+1) equidistant points coeffs = [func.subs(x, S(i) / n) for i in range(n + 1)] @@ -746,7 +721,6 @@ def tachevcoeffs(func, x, n): coeffselev = degelev(coeffselev, n // 2) return [2 * c1 - c2 for c1, c2 in zip(coeffs, coeffselev)] - def tachevcoeffsapprox(func, x, n, err=S(1) / 1000): # Same as tachevcoeffs, but coefficients are within # 'err' of the ones returned by tachevcoeffs. @@ -769,7 +743,6 @@ def tachevcoeffsapprox(func, x, n, err=S(1) / 1000): # Round to nearest multiple of delta=err/2 return [floor(v.mid() / delta + S.Half) * delta for v in ret] - def _isRandomLess(u, s): # Determines whether 'u', a uniform random variate # between 0 and 1, is less than 's'. @@ -785,7 +758,6 @@ def _isRandomLess(u, s): u[0] = (u[0] << sh) | random.randint(0, (1 << sh) - 1) u[1] = u[1] << sh - class FactoryFunction: def __init__(self, func, x): # A Bernoulli factory for continuous functions. @@ -911,7 +883,6 @@ def sim(self, coin): return 1 return 0 - class C3Function(FactoryFunction): def __init__(self, func, x, thirdderiv): # A Bernoulli factory for functions with a continuous @@ -941,7 +912,6 @@ def _errshift(self, m): def _diffwidth(self, m): return S("1.01") * 5 * self.mm * self.zz / (2 ** (2 * m + 1)) - def cheb_to_bern(n): # Transforms Chebyshev interp. coefficients to Bernstein coefficients # Rababah, Abedallah. "Transformation of Chebyshev–Bernstein polynomial basis." Computational Methods in Applied Mathematics 3.4 (2003): 608-622. @@ -958,7 +928,6 @@ def cheb_to_bern(n): ] return Matrix(mat) - def chebpoly(coeffs, x, a=0, b=1): # Polynomial on [a,b] given Chebyshev interpolant coefficients if a > b: @@ -968,7 +937,6 @@ def chebpoly(coeffs, x, a=0, b=1): for i in range(0, len(coeffs)) ) - def chebpoly2(coeffs, x): # Polynomial on [-1,1] given Chebyshev interpolant coefficients # If func^{nu} has bounded variation V, nu>=1, @@ -977,7 +945,6 @@ def chebpoly2(coeffs, x): # Theorem 7.2, Trefethen, Lloyd N., Approximation theory and approximation practice, 2013. G. Mastroianni and J. Szabados, "Jackson order of approximation by Lagrange interpolation", Acta Mathematica Hungarica, 69 (1995), 73-82. return chebpoly(coeffs, x, a=-1, b=1) - def chebdegree(eps, totvar, nu): # Upper bound on degree of polynomial to achieve error tolerance eps, # given that func^{nu} has bounded variation totvar on [-1,1]. @@ -986,14 +953,12 @@ def chebdegree(eps, totvar, nu): raise ValueError("'nu' must be 1 or greater") return ceiling(nu + (S(4 * totvar) / (pi * eps * nu)) ** (1 / S(nu))) - def chebdegree_rough(eps, totvar, nu): # Rougher upper bound on degree of polynomial to achieve error tolerance eps, # given that func^{nu} has bounded variation totvar on [-1,1]. # nu>=1 return chebdegree_01_rough(eps, totvar, nu, a=-1, b=1) - def chebdegree_01_rough(eps, totvar, nu, a=0, b=1): # Rougher upper bound on degree of polynomial to achieve error tolerance eps, # given that func^{nu} has bounded variation totvar on [a, b]. @@ -1008,7 +973,6 @@ def chebdegree_01_rough(eps, totvar, nu, a=0, b=1): ** (1 / S(nu)) ) - def chebcoeffs(func, x, n, a=0, b=1): # Chebyshev interpolant coefficients of degree n # for a continuous function func defined on [-1,1] @@ -1041,12 +1005,10 @@ def chebcoeffs(func, x, n, a=0, b=1): print(time.time() - tm) return ret - # fun=exp(-(x + S(49)/10)**2 / (S(1) / 12)) # fun1=fun.subs(x,x*24-12) # cc=chebcoeffs(fun1.subs(x,(x+1)/2),x,100) - def cheb_berncoeffs_deg(func, x, delta=S(1) / 1000, nn=10): if nn < 1: raise ValueError("'nn' must be 1 or greater") @@ -1057,7 +1019,6 @@ def cheb_berncoeffs_deg(func, x, delta=S(1) / 1000, nn=10): bern = cheb_to_bern(nn) * cc return [floor(c / delta + S.Half) * delta for c in bern] - def cheb_berncoeffs(func, x, eps=S(1) / 1000, totvar=1, nu=3): # Calculates Bernstein coefficients of a polynomial that # approximates func(x) with an error tolerance 'eps', @@ -1075,7 +1036,6 @@ def cheb_berncoeffs(func, x, eps=S(1) / 1000, totvar=1, nu=3): nn = chebdegree_01_rough(eps / 2, totvar, nu) return cheb_berncoeffs_deg(func, x, delta=delta, nn=nn) - def akahiraest(func, x, p): # Akahira et al. (1992) unbiased estimator of func(x) where # x is the expected value of a Bernoulli random variate diff --git a/moore.py b/moore.py index bd079e0b8..fb674e8af 100755 --- a/moore.py +++ b/moore.py @@ -12,7 +12,6 @@ from fractions import Fraction from interval import FInterval - class MooreSampler: """ Moore rejection sampler, for generating independent samples @@ -321,7 +320,6 @@ def _sample(self, trials=50): break return [None, trials] - if __name__ == "__main__": def bucket(v, ls, buckets): diff --git a/pcaspd.py b/pcaspd.py index 307758866..1cab1c0c4 100755 --- a/pcaspd.py +++ b/pcaspd.py @@ -704,7 +704,6 @@ [[0.287956117643289], [0.3558381900199877], [0.09809903558056805]], ] - def sRGBToSPDOtsu(srgb): """Implements Otsu and others, "Reproducing Spectral Reflectances from Tristimulus Colors", 2018.""" @@ -741,20 +740,17 @@ def sRGBToSPDOtsu(srgb): spd = SPD(refl, 10, 380) return spd - def _fileread(x): f = file(x, "r") ret = f.read() f.close() return ret - def _cmfd65(x): c = cie1931cmf(x) d = d65Illum(x) / 100.0 return [cv * d for cv in c] - def _meanarr(arrs): arrlen = len(arrs[0]) numarrs = len(arrs) @@ -764,7 +760,6 @@ def _meanarr(arrs): ret[j] += arrs[i][j] return [v * 1.0 / numarrs for v in ret] - def _pca(spdxyy, st, en): global cmf arr = [spdxyy[i][0] for i in range(st, en)] @@ -778,7 +773,6 @@ def _pca(spdxyy, st, en): cmean = matMul(cmf, mean) return [mat, minv, mean, cmean] - def _deltae(spdxyy, st, en): pcas = None try: @@ -796,7 +790,6 @@ def _deltae(spdxyy, st, en): ret += normsq return ret - def _splitAxis(spdxyy): spx = [x for x in spdxyy] spy = [x for x in spdxyy] @@ -824,7 +817,6 @@ def _splitAxis(spdxyy): return [spx[:best], spx[best:], spx[best][1][0], axis] return [spy[:best], spy[best:], spx[best][1][0], axis] - # Generates algorithm data from JSON array of reflectance curves # (stored in a file named "specs.json"). # Each curve is an array of 36 reflectance factors (from 0 through 1) diff --git a/randextract.html b/randextract.html index bdfbbe4b1..57261de79 100755 --- a/randextract.html +++ b/randextract.html @@ -118,7 +118,7 @@

      On Algorithm M

    • the chance that the first “die” shows a number less than the second “die” is the same as the chance that the first “die” shows a greater number.
    -

    More formally, P(X < Y) must be equal to P(X > Y). This relationship is equivalent to statistical indifference (Montes Gutiérrez 2014)15, (De Schuymer et al. 2003)16. This relationship works even if X and Y are dependent on each other but independent of everything else; this is easy to see if we treat X and Y as a single random “vector” [X, Y]. This is shown by the following two propositions. In the propositions below, a random variable is non-degenerate if it does not take on a single value with probability 1.

    +

    More formally, P(X < Y) must be equal to P(X > Y). This relationship is equivalent to statistical indifference (Montes Gutiérrez 2014)15, (De Schuymer et al. 2003)16. This relationship works even if X and Y are dependent on each other but independent of everything else; this is easy to see if we treat X and Y as a single random “vector” [X, Y]. This is shown by the following two propositions. In the propositions below, a random variable is nondegenerate if it does not take on a single value with probability 1.

    Proposition 1. Let X and Y be real-valued nondegenerate random variables. Then Algorithm M outputs 1 or 0 with equal probability if and only if X and Y are statistically indifferent.

    diff --git a/randomfunc.html b/randomfunc.html index 11a710c04..fcb58fe0c 100755 --- a/randomfunc.html +++ b/randomfunc.html @@ -1013,9 +1013,9 @@

    A Note on Sorting Random Variates

    -

    General Non-Uniform Distributions

    +

    General Non-Uniform Distributions

    -

    Some applications need to choose random values such that some of them have a greater chance to be chosen than others (a non-uniform distribution). Most of the techniques in this section show how to use the uniform random integer methods to generate such random values.

    +

    Some applications need to choose random values such that some of them have a greater chance to be chosen than others (a nonuniform distribution). Most of the techniques in this section show how to use the uniform random integer methods to generate such random values.

    @@ -1325,7 +1325,7 @@

    Transformations of Random Variates

    -

    Specific Non-Uniform Distributions

    +

    Specific Non-Uniform Distributions

    This section contains information on some of the most common nonuniform probability distributions.

    @@ -2270,7 +2270,7 @@

    Specific Distributions

    -

    Index of Non-Uniform Distributions

    +

    Index of Non-Uniform Distributions

    Many distributions here require random real numbers.

    @@ -2557,7 +2557,7 @@

    Random Points Inside a Bo
    1. The Python sample code contains a method for generating a random point on the surface of an ellipsoid modeling the Earth.
    2. Sampling a half-ball, half-sphere, or half-shell can be done by sampling a full ball or shell and replacing one of the dimensions of the result with its absolute value.
    3. -
    4. Lacko and Harman (2012)111 defined a family of non-uniform distributions of points inside a ball: generate RandomPointOnSphere(dims, r*pow(BetaDist(dims/p, d/p), 1.0/p),p) where r>0 is the radius, dims and p are as in RandomPointOnSphere, and d>=0 is a shape parameter. If d = p, the distribution is uniform in the ball.
    5. +
    6. Lacko and Harman (2012)111 defined a family of nonuniform distributions of points inside a ball: generate RandomPointOnSphere(dims, r*pow(BetaDist(dims/p, d/p), 1.0/p),p) where r>0 is the radius, dims and p are as in RandomPointOnSphere, and d>=0 is a shape parameter. If d = p, the distribution is uniform in the ball.
    diff --git a/randomgendoc.md b/randomgendoc.md index 76437a594..f7f75ed3e 100755 --- a/randomgendoc.md +++ b/randomgendoc.md @@ -1279,7 +1279,7 @@ CLASSES | wiener(self, st, en, step=1.0, mu=0.0, sigma=1.0) | Generates random variates following a Wiener | process (Brownian motion). Each element of the return - | value contains a timestamp and a random variate in that order. + | value contains a time stamp and a random variate in that order. | | zero_or_one(self, px, py) | Returns 1 at probability px/py, 0 otherwise. @@ -1718,7 +1718,7 @@ CLASSES | - Goyal, V. and Sigman, K. 2012. On simulating a class of Bernstein | polynomials. ACM Transactions on Modeling and Computer Simulation 22(2), | Article 12 (March 2012), 5 pages. - | - Giulio Morina. Krzysztof Łatuszyński. Piotr Nayar. Alex Wendland. "From the Bernoulli factory to a dice enterprise via perfect sampling of Markov chains." Ann. Appl. Probab. 32 (1) 327 - 359, February 2022. https://doi.org/10.1214/21-AAP1679 + | - Giulio Morina. Krzysztof Łatuszyński. Piotr Nayar. Alex Wendland. "From the Bernoulli factory to a dice enterprise via perfect sampling of Markov chains." Ann. Appl. Probab. 32 (1) 327 - 359, February 2022. [https://doi.org/10.1214/21-AAP1679](https://doi.org/10.1214/21-AAP1679) | - Dughmi, Shaddin, Jason Hartline, Robert D. Kleinberg, and Rad Niazadeh. "Bernoulli factories and black-box reductions in mechanism design." Journal of the ACM (JACM) 68, no. 2 (2021): 1-30. | - Gonçalves, F. B., Łatuszyński, K. G., Roberts, G. O. (2017). Exact Monte | Carlo likelihood-based inference for jump-diffusion processes. @@ -4491,4 +4491,3 @@ FILE /home/peter/Documents/SharpDevelopProjects/peteroupc.github.io/betadist.py ``` - diff --git a/sitemap.xml b/sitemap.xml index 9ccd73acd..52f6d0c27 100755 --- a/sitemap.xml +++ b/sitemap.xml @@ -137,7 +137,7 @@ https://peteroupc.github.io/MailLib/api/com.upokecenter.mail.DataUris.html2024-04-27T05:23:05Z https://peteroupc.github.io/MailLib/api/com.upokecenter.text.ICharacterEncoding.html2021-02-15T12:41:19Z https://peteroupc.github.io/MailLib/api/com.upokecenter.mail.LanguageTags.StringAndQuality.html2021-02-15T12:41:19Z -https://peteroupc.github.io/MailLib/api/com.upokecenter.mail.MediaType.html2024-12-26T01:03:59Z +https://peteroupc.github.io/MailLib/api/com.upokecenter.mail.MediaType.html2024-12-26T09:03:57Z https://peteroupc.github.io/MailLib/api/com.upokecenter.mail.Address.html2024-04-27T05:23:05Z https://peteroupc.github.io/MailLib/api/com.upokecenter.text.Encodings.html2021-02-15T12:41:19Z https://peteroupc.github.io/MailLib/api/com.upokecenter.mail.DataUrls.html2024-04-27T05:23:05Z @@ -148,7 +148,7 @@ https://peteroupc.github.io/MailLib/api/com.upokecenter.mail.transforms.QuotedPrintableTransform.html2021-02-15T12:41:19Z https://peteroupc.github.io/MailLib/api/com.upokecenter.mail.transforms.Base64Transform.html2021-02-15T12:41:19Z https://peteroupc.github.io/MailLib/api/com.upokecenter.mail.MediaType.Builder.html2024-12-24T16:28:42Z -https://peteroupc.github.io/MailLib/api/com.upokecenter.text.NormalizerInput.html2024-04-27T05:23:05Z +https://peteroupc.github.io/MailLib/api/com.upokecenter.text.NormalizerInput.html2024-12-26T09:03:57Z https://peteroupc.github.io/MailLib/api/com.upokecenter.mail.transforms.BoundaryCheckerTransform.html2021-02-15T12:41:19Z https://peteroupc.github.io/MailLib/api/com.upokecenter.text.ICharacterInput.html2021-02-15T12:41:19Z https://peteroupc.github.io/MailLib/api/com.upokecenter.mail.transforms.LiberalSevenBitTransform.html2021-02-15T12:41:19Z @@ -159,7 +159,7 @@ https://peteroupc.github.io/requests.html2024-12-24T16:29:01Z https://peteroupc.github.io/exporand.pdf2024-11-09T22:04:16Z https://peteroupc.github.io/logconcave.py2024-12-24T11:58:23Z -https://peteroupc.github.io/colorutil.zip2024-12-26T01:05:09Z +https://peteroupc.github.io/colorutil.zip2024-12-26T09:04:31Z https://peteroupc.github.io/TurtleParser/index.html2024-08-08T10:19:24Z https://peteroupc.github.io/TurtleParser/api/com.upokecenter.util.IWriter.html2021-02-15T12:41:19Z https://peteroupc.github.io/TurtleParser/api/com.upokecenter.util.ArrayWriter.html2021-02-15T12:41:19Z @@ -483,7 +483,7 @@ https://peteroupc.github.io/html3dutil/H3DU.ArcLengthParamCurve.html2024-12-23T00:49:50Z https://peteroupc.github.io/html3dutil/H3DU.LightSource.html2024-12-22T08:56:26Z https://peteroupc.github.io/html3dutil/H3DU.Semantic.html2024-12-07T16:29:50Z -https://peteroupc.github.io/html3dutil/extras_contourlines.html2024-12-26T01:04:05Z +https://peteroupc.github.io/html3dutil/extras_contourlines.html2024-12-26T09:04:00Z https://peteroupc.github.io/html3dutil/H3DU.FrameBufferInfo.html2023-03-20T06:26:05Z https://peteroupc.github.io/html3dutil/extras_stl.html2024-12-25T12:42:13Z https://peteroupc.github.io/html3dutil/H3DU.CurveBuilder.html2024-12-25T12:42:09Z @@ -539,7 +539,7 @@ https://peteroupc.github.io/Numbers/docs/PeterO.Numbers.EDecimals.html2024-12-24T11:58:53Z https://peteroupc.github.io/Numbers/api/index.html2024-04-27T05:23:06Z https://peteroupc.github.io/Numbers/api/com.upokecenter.numbers.EDecimal.html2024-12-24T11:58:04Z -https://peteroupc.github.io/Numbers/api/com.upokecenter.numbers.EFloat.html2024-12-26T01:04:01Z +https://peteroupc.github.io/Numbers/api/com.upokecenter.numbers.EFloat.html2024-12-26T09:03:58Z https://peteroupc.github.io/Numbers/api/Home.html2024-04-27T05:23:06Z https://peteroupc.github.io/Numbers/api/com.upokecenter.numbers.ETrapException.html2024-04-27T05:23:06Z https://peteroupc.github.io/Numbers/api/com.upokecenter.numbers.EContext.html2024-12-24T11:58:53Z @@ -550,13 +550,13 @@ https://peteroupc.github.io/Numbers/api/com.upokecenter.numbers.EDecimalExtras.html2021-02-15T12:41:19Z https://peteroupc.github.io/Numbers/api/com.upokecenter.numbers.EDecimals.html2024-12-24T16:28:43Z https://peteroupc.github.io/Numbers/api/com.upokecenter.numbers.ERational.html2024-12-24T11:58:05Z -https://peteroupc.github.io/exporand.html2024-12-26T01:04:35Z +https://peteroupc.github.io/exporand.html2024-12-26T09:04:15Z https://peteroupc.github.io/randomgen.py2024-12-24T11:58:30Z https://peteroupc.github.io/bernoullicorrect.pdf2024-11-09T22:03:42Z https://peteroupc.github.io/hash.pdf2024-11-09T22:04:26Z https://peteroupc.github.io/uniformsum.html2024-04-27T05:23:14Z https://peteroupc.github.io/Enriched/index.html2023-03-20T06:26:05Z -https://peteroupc.github.io/randomgen.zip2024-12-26T01:05:09Z +https://peteroupc.github.io/randomgen.zip2024-12-26T09:04:31Z https://peteroupc.github.io/jump.html2024-12-24T11:58:55Z https://peteroupc.github.io/morealg.html2024-04-27T05:23:14Z https://peteroupc.github.io/randomtest.html2024-12-24T11:58:55Z @@ -571,7 +571,7 @@ https://peteroupc.github.io/requestsother.html2024-12-24T16:29:01Z https://peteroupc.github.io/requestsother.pdf2024-11-09T22:05:30Z https://peteroupc.github.io/Fuzzer/index.html2024-12-22T08:56:27Z -https://peteroupc.github.io/randomfunc.html2024-12-26T01:04:33Z +https://peteroupc.github.io/randomfunc.html2024-12-26T09:04:14Z https://peteroupc.github.io/hqprng.html2024-12-24T11:58:56Z https://peteroupc.github.io/filenames.pdf2024-11-09T22:04:19Z https://peteroupc.github.io/colorpicker/index.html2024-12-23T00:49:51Z @@ -579,9 +579,9 @@ https://peteroupc.github.io/jump.pdf2024-11-09T22:04:36Z https://peteroupc.github.io/svg.html2024-12-23T00:49:53Z https://peteroupc.github.io/graphics.pdf2024-11-09T22:04:23Z -https://peteroupc.github.io/randextract.html2024-12-26T01:04:36Z +https://peteroupc.github.io/randextract.html2024-12-26T09:04:15Z https://peteroupc.github.io/svg.pdf2024-11-09T22:05:41Z -https://peteroupc.github.io/bernoulli.zip2024-12-26T01:05:10Z +https://peteroupc.github.io/bernoulli.zip2024-12-26T09:04:32Z https://peteroupc.github.io/petero-csharp/index.html2024-12-22T08:56:27Z https://peteroupc.github.io/petero-csharp/docs/index.html2024-04-27T05:23:10Z https://peteroupc.github.io/petero-csharp/docs/PeterO.InvariantText.html2024-04-27T05:23:10Z @@ -606,5 +606,5 @@ https://peteroupc.github.io/petero-csharp/docs/PeterO.IniEntry.html2024-04-27T05:23:10Z https://peteroupc.github.io/petero-csharp/docs/PeterO.ArrayUtil.html2024-04-27T05:23:10Z https://peteroupc.github.io/colorgen.pdf2024-11-09T22:04:03Z -https://peteroupc.github.io/betadist.py2024-12-26T01:04:58Z +https://peteroupc.github.io/betadist.py2024-12-26T09:04:25Z