Ticket #29230: sigma_gauss_function.py

File sigma_gauss_function.py, 3.0 KB (added by gh-garambois, 3 years ago)
Line 
1r"""
2sigma_gauss function (sum of divisors) applicable to a Gaussian integer.
3
4This function is an extension of the sigma function in number theory to Gaussian integers.
5
6
7INPUT:
8
9    ''z'' -- GaussianInteger
10
11OUTPUT:
12
13    ''sigma_gauss(z)'' -- GaussianInteger, the sum of the Gaussian integer divisors of z
14
15    Caution: The Gaussian integer prime factors are only taken in the first quadrant,
16    see why in the references below.
17
18EXAMPLES:
19
20The sum of the divisors of Gauss integers z = 2*I + 2, z = 3*I + 2, z = 13::
21
22    sage: sigma_gauss(2*I + 2)
23    5*I
24    sage: sigma_gauss(3*I + 2)
25    3*I + 3
26    sage: sigma_gauss(5)
27    8*I + 4
28
29AUTHORS:
30
31    Paul Zimmermann and Jean-Luc Garambois (2019): initial version
32
33REFERENCES:
34
35    [1] http://mathworld.wolfram.com/DivisorFunction.html
36        Equivalent on the Mathematica software to DivisorSigma [1, z, GaussianIntegers -> True]
37    [2] https://www.jstor.org/stable/2312472?seq=1
38    [3] https://encompass.eku.edu/etd/158/
39"""
40# ****************************************************************************
41#       Copyright (C) 2019 Jean-Luc Garambois <jlgarambois@gmail.com>
42#
43# This program is free software: you can redistribute it and/or modify
44# it under the terms of the GNU General Public License as published by
45# the Free Software Foundation, either version 2 of the License, or
46# (at your option) any later version.
47#                  https://www.gnu.org/licenses/
48# ****************************************************************************
49
50from sage.functions.other import real, imag
51from sage.rings.number_field.order import GaussianIntegers
52
53R = GaussianIntegers()
54
55def sigma_gauss_aux (l):
56    s = 1
57    e = 0
58    while e < len(l):
59        p = l[e][0]
60        k = l[e][1]
61        s = s * (p**(k+1) - 1) / (p-1)
62        e = e + 1
63    return R(s)
64
65def sigma_gauss(z):
66    ll = list(R(z).factor())
67    l = []
68    i = 0
69
70# ****************************************************************************
71# Gaussian prime factors must be selected in the first quadrant.
72# ****************************************************************************
73
74    while i < len(ll):
75        nz = ll[i][0]
76        if real(nz) < 0 and imag(nz) < 0:
77            nz = -nz
78        if real(nz) < 0 and imag(nz) >= 0:
79            nz = -I * nz
80        if real(nz) >= 0 and imag(nz) < 0:
81            nz = I * nz
82        assert real(nz) >= 0 and imag(nz) >= 0
83        lz = [nz, ll[i][1]]
84        l.append(lz)
85        i = i + 1
86    r = 0
87
88# ****************************************************************************
89# Be careful not to have the same factor twice after their selection
90# in the first quadrant and before applying the "sigma_gauss_aux" function.
91# ****************************************************************************
92
93    while r < len(l) - 1:
94        t = r + 1
95        while t < len(l):
96            if l[r][0] == l[t][0]:
97                l[r][1] = l[r][1] + l[t][1]
98                del l[t]
99                t = t - 1
100            t = t + 1
101        r= r + 1
102    for x, e in l:
103        assert e != 0
104    return sigma_gauss_aux(l)