Opened 10 years ago
Last modified 8 years ago
#12418 closed enhancement
adding Delsarte bound for codes — at Version 9
Reported by: | dimpase | Owned by: | wdj |
---|---|---|---|
Priority: | major | Milestone: | sage-5.12 |
Component: | coding theory | Keywords: | |
Cc: | jpang, kini, wdj, ncohen, ppurka, ptrrsn_1 | Merged in: | |
Authors: | Reviewers: | ||
Report Upstream: | N/A | Work issues: | |
Branch: | Commit: | ||
Dependencies: | Stopgaps: |
Description (last modified by )
Delsarte bound for codes, aka Linear Programming bound, is easy to implement in Sage. Here is a quick and dirty code that does it; so the problem would be to integrate this properly.
def Kra(n,q,l,i): # K^{n,q}_l(i), i.e Krawtchouk polynomial return sum([((-1)**j)*((q-1)**(l-j))*binomial(i,j)*binomial(n-i,l-j) for j in range(l+1)]) def roundres(x): # this is a quick and unsafe way to round the result... import math tol = 0.0001 if math.ceil(x)-x < tol: return int(math.ceil(x)) if x-math.floor(x) < tol: return int(math.floor(x)) return x # @cached_function def delsarte_bound(n, q, d, d_star=1, q_base=0, return_log=True,\ isinteger=False, return_data=False): p = MixedIntegerLinearProgram(maximization=True) A = p.new_variable(integer=isinteger) # A>=0 is assumed p.set_objective(sum([A[r] for r in range(n+1)])) p.add_constraint(A[0]==1) for i in range(1,d): p.add_constraint(A[i]==0) for j in range(1,n+1): rhs = sum([Kra(n,q,j,r)*A[r] for r in range(n+1)]) if j >= d_star: p.add_constraint(0*A[0] <= rhs) else: # rhs is proportional to j-th weight of the dual code p.add_constraint(0*A[0] == rhs) try: bd=p.solve() except sage.numerical.mip.MIPSolverException, exc: print "Solver exception: ", exc, exc.args if return_data: return A,p,False return False if q_base > 0: # rounding the bound down to the nearest power of q_base, # for q=q_base^m # qb = factor(q).radical() # if len(qb) == 1: # base = qb.expand() # bd = base^int(log(bd, base=base)) if return_log: # bd = int(log(roundres(bd), base=q_base)) # unsafe: # loss of precision bd = roundres(log(bd, base=q_base)) else: # bd = q_base^int(log(roundres(bd), base=q_base)) bd = q_base^roundres(log(bd, base=q_base)) if return_data: return A,p,bd return bd
d_star is the minimal distance of the dual code (if it exists at all)
If q_base is 0, just compute the upper bound.
If q_base is >0, it is assumed to be a prime power, and the code assumed
to be additive over this field (i.e. the dual code exists,and its weight enumerator is
obtained by applying MacWilliams transform
--- the matrix A of the LP times the
weight enumerator of our code), then the output is the corr. dimension, i.e.
floor(log(bound, q_base))
.
One obstacle for this to work well in big dimensions is a lack of arbitrary precision LP solver backend available in Sage. This is (almost - i.e. the corresponding ticket is still not 100% ready, as reviewers think) taken care of by #12533, which should be made a dependence for this ticket.
Change History (9)
comment:1 Changed 10 years ago by
comment:2 Changed 10 years ago by
- Cc kini added
comment:3 Changed 10 years ago by
- Cc wdj added
comment:4 Changed 10 years ago by
- Cc ncohen added
- Description modified (diff)
comment:5 Changed 10 years ago by
comment:6 Changed 10 years ago by
- Cc ppurka added
comment:7 Changed 10 years ago by
The function named delsarte_bound
should be renamed to something like delsarte_bound_hamming_space
. This is so that in future other functions like delsarte_bound_johnson_space
, delsarte_bound_permutation_space
, etc can be added easily, without having inconsistencies in naming.
comment:8 Changed 10 years ago by
- Cc ptrrsn_1 added
comment:9 Changed 10 years ago by
- Description modified (diff)
(GLPK can solve non-integer rational LP. It is not exposed, but may not be too hard either)