2
2
Optimized counting of congruence solutions
3
3
"""
4
4
from sage.arith.misc import is_prime, kronecker as kronecker_symbol, valuation
5
- from sage.rings.finite_rings.integer_mod cimport IntegerMod_gmp
6
5
from sage.rings.finite_rings.integer_mod import Mod
7
6
from sage.rings.finite_rings.integer_mod_ring import IntegerModRing
8
7
@@ -103,11 +102,6 @@ cdef CountAllLocalTypesNaive_cdef(Q, p, k, m, zvec, nzvec):
103
102
R = p ** k
104
103
Q1 = Q.change_ring(IntegerModRing(R))
105
104
106
- # Cython Variables
107
- cdef IntegerMod_gmp zero, one
108
- zero = IntegerMod_gmp(IntegerModRing(R), 0 )
109
- one = IntegerMod_gmp(IntegerModRing(R), 1 )
110
-
111
105
# Initialize the counting vector
112
106
count_vector = [0 for i in range (6 )]
113
107
@@ -119,7 +113,7 @@ cdef CountAllLocalTypesNaive_cdef(Q, p, k, m, zvec, nzvec):
119
113
m1 = Mod(m, R)
120
114
121
115
# Count the local solutions
122
- for i from 0 <= i < R_n:
116
+ for i in range ( R_n) :
123
117
124
118
# Perform a carry (when value = R-1) until we can increment freely
125
119
ptr = len (v)
@@ -128,20 +122,20 @@ cdef CountAllLocalTypesNaive_cdef(Q, p, k, m, zvec, nzvec):
128
122
ptr += - 1
129
123
130
124
# Only increment if we're not already at the zero vector =)
131
- if ( ptr > 0 ) :
125
+ if ptr > 0 :
132
126
v[ptr- 1 ] += 1
133
127
134
128
# Evaluate Q(v) quickly
135
129
tmp_val = Mod(0 , R)
136
130
for a from 0 <= a < n:
137
131
for b from a <= b < n:
138
- tmp_val += Q1[a,b] * v[a] * v[b]
132
+ tmp_val += Q1[a, b] * v[a] * v[b]
139
133
140
134
# Sort the solution by it's type
141
- # if ( Q1(v) == m1) :
142
- if ( tmp_val == m1) :
135
+ # if Q1(v) == m1:
136
+ if tmp_val == m1:
143
137
solntype = local_solution_type_cdef(Q1, p, v, zvec, nzvec)
144
- if ( solntype != 0 ) :
138
+ if solntype != 0 :
145
139
count_vector[solntype] += 1
146
140
147
141
# Generate the Bad-type and Total counts
@@ -187,11 +181,6 @@ def CountAllLocalTypesNaive(Q, p, k, m, zvec, nzvec):
187
181
return CountAllLocalTypesNaive_cdef(Q, p, k, m, zvec, nzvec)
188
182
189
183
190
-
191
-
192
-
193
-
194
-
195
184
cdef local_solution_type_cdef(Q, p, w, zvec, nzvec):
196
185
"""
197
186
Internal routine to check if a given solution vector `w` (of `Q(w) =
@@ -209,91 +198,87 @@ cdef local_solution_type_cdef(Q, p, w, zvec, nzvec):
209
198
210
199
# Check if the solution satisfies the zvec "zero" congruence conditions
211
200
# (either zvec is empty or its components index the zero vector mod p)
212
- if ( zvec is None ) or ( not zvec) :
201
+ if zvec is None or not zvec:
213
202
zero_flag = True
214
203
else :
215
204
zero_flag = False
216
205
i = 0
217
- while ( ( i < len (zvec)) and (( w[zvec[i]] % p) == 0 ) ) : # Increment so long as our entry is zero (mod p)
206
+ while i < len (zvec) and not w[zvec[i]] % p: # Increment so long as our entry is zero (mod p)
218
207
i += 1
219
- if ( i == len (zvec) ): # If we make it through all entries then the solution is zero (mod p)
208
+ if i == len (zvec): # If we make it through all entries then the solution is zero (mod p)
220
209
zero_flag = True
221
210
222
-
223
211
# DIAGNOSTIC
224
- # print("IsLocalSolutionType: Finished the Zero congruence condition test \n")
212
+ # print("IsLocalSolutionType: Finished the Zero congruence condition test \n")
225
213
226
214
if not zero_flag:
227
215
return < long > 0
228
216
229
217
# DIAGNOSTIC
230
- # print("IsLocalSolutionType: Passed the Zero congruence condition test \n")
218
+ # print("IsLocalSolutionType: Passed the Zero congruence condition test \n")
231
219
232
220
# Check if the solution satisfies the nzvec "nonzero" congruence conditions
233
221
# (nzvec is non-empty and its components index a non-zero vector mod p)
234
- if ( nzvec is None ) :
222
+ if nzvec is None :
235
223
nonzero_flag = True
236
- elif ( len (nzvec) == 0 ) :
224
+ elif len (nzvec) == 0 :
237
225
nonzero_flag = False # Trivially no solutions in this case!
238
226
else :
239
227
nonzero_flag = False
240
228
i = 0
241
- while (( not nonzero_flag) and ( i < len (nzvec)) ):
242
- if (( w[nzvec[i]] % p) ! = 0 ) :
229
+ while not nonzero_flag and i < len (nzvec):
230
+ if w[nzvec[i]] % p:
243
231
nonzero_flag = True # The non-zero condition is satisfied when we find one non-zero entry
244
232
i += 1
245
233
246
234
if not nonzero_flag:
247
235
return < long > 0
248
236
249
-
250
237
# Check if the solution has the appropriate (local) type:
251
238
# -------------------------------------------------------
252
239
253
240
# 1: Check Good-type
254
- for i from 0 <= i < n :
255
- if ((( w[i] % p) ! = 0 ) and (( Q[i,i] % p) ! = 0 )) :
241
+ for i in range (n) :
242
+ if w[i] % p and Q[i, i] % p:
256
243
return < long > 1
257
- if ( p == 2 ) :
258
- for i from 0 <= i < (n - 1 ):
259
- if ((( Q[i,i+ 1 ] % p) ! = 0 ) and (((w[i] % p) ! = 0 ) or ((w[i+ 1 ] % p) ! = 0 ))):
244
+ if p == 2 :
245
+ for i in range (n - 1 ):
246
+ if ((Q[i, i+ 1 ] % p)) and (((w[i] % p)) or ((w[i+ 1 ] % p))):
260
247
return < long > 1
261
248
262
-
263
249
# 2: Check Zero-type
264
250
Zero_flag = True
265
- for i from 0 <= i < n :
266
- if (( w[i] % p) ! = 0 ) :
251
+ for i in range (n) :
252
+ if w[i] % p:
267
253
Zero_flag = False
268
254
if Zero_flag:
269
255
return < long > 2
270
256
271
-
272
257
# Check if wS1 is zero or not
273
258
wS1_nonzero_flag = False
274
259
for i from 0 <= i < n:
275
260
276
261
# Compute the valuation of each index, allowing for off-diagonal terms
277
- if ( Q[i,i] == 0 ) :
278
- if ( i == 0 ) :
279
- val = valuation(Q[i,i+ 1 ], p) # Look at the term to the right
280
- elif ( i == n - 1 ) :
281
- val = valuation(Q[i- 1 ,i], p) # Look at the term above
262
+ if Q[i, i] == 0 :
263
+ if i == 0 :
264
+ val = valuation(Q[i, i+ 1 ], p) # Look at the term to the right
265
+ elif i == n - 1 :
266
+ val = valuation(Q[i- 1 , i], p) # Look at the term above
282
267
else :
283
- val = valuation(Q[i,i+ 1 ] + Q[i- 1 ,i], p) # Finds the valuation of the off-diagonal term since only one isn't zero
268
+ val = valuation(Q[i, i+ 1 ] + Q[i- 1 , i], p) # Finds the valuation of the off-diagonal term since only one isn't zero
284
269
else :
285
- val = valuation(Q[i,i], p)
270
+ val = valuation(Q[i, i], p)
286
271
287
272
# Test each index
288
- if (( val == 1 ) and (( w[i] % p) ! = 0 )) :
273
+ if val == 1 and w[i] % p:
289
274
wS1_nonzero_flag = True
290
275
291
276
# 4: Check Bad-type I
292
- if ( wS1_nonzero_flag is True ) :
277
+ if wS1_nonzero_flag:
293
278
return < long > 4
294
279
295
280
# 5: Check Bad-type II
296
- if ( wS1_nonzero_flag is False ) :
281
+ if wS1_nonzero_flag:
297
282
return < long > 5
298
283
299
284
# Error if we get here! =o
0 commit comments