Ticket #29230: sigma_gauss_function.sage

File sigma_gauss_function.sage, 1.4 KB (added by gh-garambois, 3 years ago)
Line 
1#Calculation program sigma(z), z Gaussian integer : calculates the sum of the divisors of z with the Gaussian integer prime factors taken in the first quadrant.
2#References: equivalent to the function DivisorSigma [1, z, GaussianIntegers -> True] in the Mathematica software: http://mathworld.wolfram.com/DivisorFunction.html
3#References: https://www.jstor.org/stable/2312472?seq=1
4#References: https://encompass.eku.edu/etd/158/
5
6R = GaussianIntegers()
7
8def sigma_gauss_aux (l):
9   s = 1
10   e=0
11   while e<len(l):
12      p = l[e][0]
13      k = l[e][1]
14      s = s * (p^(k+1)-1) / (p-1)
15      e+=1
16   return R(s)
17
18def sigma_gauss(z):
19   ll=list(R(z).factor())
20   l=[]
21   i=0
22
23#Gaussian prime factors must be selected in the first quadrant.
24
25   while i<len(ll):
26      nz=ll[i][0]
27      if real(nz)<0 and imag(nz)<0:
28         nz=-nz
29      if real(nz)<0 and imag(nz)>=0:
30         nz=-I*nz
31      if real(nz)>=0 and imag(nz)<0:
32         nz=I*nz
33      assert real(nz)>=0 and imag(nz)>=0
34      lz=[nz,ll[i][1]]
35      l.append(lz)
36      i+=1
37   r=0
38
39#Be careful not to have the same factor twice after their selection in the first quadrant and before applying the "sigma_gauss_aux" formula.
40
41   while r<len(l)-1:
42      t=r+1
43      while t<len(l):
44         if l[r][0]==l[t][0]:
45            l[r][1]=l[r][1]+l[t][1]
46            del l[t]
47            t-=1
48         t+=1
49      r+=1
50   for x,e in l:
51      assert e != 0
52   return sigma_gauss_aux(l)