Ticket #4333: trac_4333_remove_bernoulli_py.patch

File trac_4333_remove_bernoulli_py.patch, 4.9 KB (added by mhampton, 5 years ago)

based against 3.2 alpha0

  • sage/rings/arith.py 

    # HG changeset patch
# User Marshall Hampton <hamptonio@gmail.com>
# Date 1224701689 18000
# Node ID a33bab2bb3f5e22875ed15c70f685f0d5b27618e
# Parent  a7b1012041c9ac599cc11ba38522d40c58701362
removed bernoulli.py since it was broken and redundant

diff -r a7b1012041c9 -r a33bab2bb3f5 sage/rings/arith.py
    a b  
    161161            'gap'  -- use GAP 
    162162            'gp'   -- use PARI/GP interpreter 
    163163            'magma' -- use MAGMA (optional) 
    164             'python' -- use pure Python implementation 
    165164            'bernmm' -- use bernmm package (a multimodular algorithm) 
    166165        num_threads -- positive integer, number of threads to use 
    167166                       (only used for bernmm algorithm) 
    … …  
    181180        -691/2730 
    182181        sage: bernoulli(12, algorithm='pari')     
    183182        -691/2730 
    184         sage: bernoulli(12, algorithm='python') 
    185         -691/2730 
    186183        sage: bernoulli(12, algorithm='bernmm') 
    187184        -691/2730 
    188185        sage: bernoulli(12, algorithm='bernmm', num_threads=4) 
    189186        -691/2730 
    190187 
     188    TESTS: 
     189        sage: algs = ['gap','gp','pari','bernmm'] #long time 
     190        sage: vals = [[bernoulli(i,algorithm = j) for j in algs] for i in range(2,2255)] #long time 
     191        sage: union([len(union(x))==1 for x in vals]) #long time 
     192        [True] 
     193        sage: algs = ['gp','pari','bernmm'] #long time 
     194        sage: vals = [[bernoulli(i,algorithm = j) for j in algs] for i in range(2256,5000)] #long time 
     195        sage: union([len(union(x))==1 for x in vals]) #long time 
     196        [True] 
    191197    \note{If $n>50000$ then algorithm = 'gp' is used instead of 
    192198    algorithm = 'pari', since the C-library interface to PARI 
    193199    is limited in memory for individual operations.} 
    … …  
    218224        import sage.interfaces.gp         
    219225        x = sage.interfaces.gp.gp('bernfrac(%s)'%n) 
    220226        return Rational(x) 
    221     elif algorithm == 'python': 
    222         import sage.rings.bernoulli 
    223         return sage.rings.bernoulli.bernoulli_python(n) 
    224227    elif algorithm == 'bernmm': 
    225228        import sage.rings.bernmm 
    226229        return sage.rings.bernmm.bernmm_bern_rat(n, num_threads) 

  • deleted file sage/rings/bernoulli.py 

    diff -r a7b1012041c9 -r a33bab2bb3f5 sage/rings/bernoulli.py
    + -  
    1 """ 
    2 Computation of Bernoulli numbers 
    3 """ 
    4  
    5 ########################################################################### 
    6 #       Copyright (C) 2006 William Stein <wstein@gmail.com> 
    7 # 
    8 #  Distributed under the terms of the GNU General Public License (GPL) 
    9 #                  http://www.gnu.org/licenses/ 
    10 ########################################################################### 
    11  
    12  
    13 def bernoulli_python(m): 
    14     r""" 
    15     Returns the Bernoulli number $B_m$ computed using 
    16     a Python implementation. 
    17     """ 
    18     import sage.rings.all as rings 
    19     import sage.misc.all as misc 
    20  
    21     m = int(m) 
    22     if m == 1: 
    23         return rings.Rational('-1/2') 
    24     if m % 2: 
    25         return 0 
    26      
    27     tm = misc.verbose('bernoulli_python setup') 
    28     RR = rings.RealField(4*4*m) 
    29     tm = misc.verbose('computing pi...', tm) 
    30     pi = RR.pi() 
    31     tm = misc.verbose('computing factorial...', tm) 
    32     m_factorial = rings.factorial(m) 
    33     tm = misc.verbose('computing pi pow...', tm) 
    34     K = 2*m_factorial/((2*pi)**m) 
    35     tm = misc.verbose('computing P...', tm) 
    36     P = rings.prime_range(m+2) 
    37  
    38     # IDIOTIC -- should compute by factoring m! 
    39     d = misc.prod([p for p in P if m % (p-1) == 0]) 
    40      
    41     tm = misc.verbose('computing N...', tm) 
    42     #N = long((K*d)**(1.0/(m-1.0)) + 1) 
    43     dK = d*K 
    44     n = str(dK).find('.') + 1 
    45     N = 10**(float(n)/float(m-1)) 
    46     tm = misc.verbose('N = %s; computing product...'%N, tm) 
    47     assert N < m 
    48     z = 1 
    49     z_prev = z 
    50     Rl = rings.RealField() 
    51     for p in P: 
    52         if p > N: 
    53             break 
    54         z /= 1 - 1/(RR(p)**RR(m)) 
    55         diff = abs(z - z_prev).log() 
    56         if diff < -16*m: 
    57             break 
    58         z_prev =z  
    59     a = long((dK*z).ceil()) 
    60     if m % 4 == 0: 
    61         a = -a 
    62     return rings.Rational(a)/rings.Rational(d) 
    63  
    64  
    65  
    66  
    67 def bernoulli_cf(m): 
    68     r""" 
    69     Returns the Bernoulli number $B_m$. 
    70     """ 
    71     import sage.rings.all as rings 
    72     import sage.misc.all as misc 
    73  
    74     m = int(m) 
    75     if m == 1: 
    76         return rings.Rational('-1/2') 
    77     if m % 2: 
    78         return 0 
    79      
    80     tm = misc.verbose('bernoulli_python setup') 
    81     RR = rings.RealField(4*4*m) 
    82     tm = misc.verbose('computing pi...', tm) 
    83     pi = RR.pi() 
    84     tm = misc.verbose('computing factorial...', tm) 
    85     m_factorial = rings.factorial(m) 
    86     tm = misc.verbose('computing pi pow...', tm) 
    87     pi_pow = (2*pi)**m 
    88     tm = misc.verbose('computing quo...', tm) 
    89     K = 2*m_factorial/pi_pow 
    90     tm = misc.verbose('computing zeta times K...', tm) 
    91     s = RR(m).zeta() * K 
    92     tm = misc.verbose('done', tm) 
    93     return s