# Copyright (c) IMToolkit Development Team
# This toolkit is released under the MIT License, see LICENSE.txt
import os
import sys
import glob
import re
import time
import shutil
from scipy import special
import numpy as np
import itertools
from imtoolkit import *
[docs]def getHammingDistanceTable(MCK, indsdec):
# This method is not dependent on Q
hds = np.zeros((MCK, MCK), dtype = np.int)
imax = MCK * (MCK - 1) / 2
i = 0
for y in range(MCK):
for x in range(y + 1, MCK):
hammingdis = countErrorBits(indsdec[y], indsdec[x])
hds[y][x] = hds[x][y] = hammingdis
i += 1
print("%.3f percent completed." % (i / imax * 100.0))
return hds
[docs]def getGoodDecsTable(M, K):
inds = list(itertools.combinations(range(M), K))
indsv = convertIndsToVector(inds, M)
indsdec = convertIndsToIndsDec(inds, M)
MCK = len(inds)
hds = np.zeros((MCK, MCK), dtype = np.int) # This is not dependent on Q
minHT = 4
newdecs = {}
for x in range(1, MCK):
hds[0][x] = np.sum(np.logical_xor(indsv[0].reshape(M), indsv[x].reshape(M)))
#print(hds[0])
while True:
gooddecs = where(hds[0] >= minHT)[0].tolist()
if len(gooddecs) == 0:
break
print("minHT = %d" % (minHT))
newdecs[minHT] = [0]
newdecs[minHT].extend(gooddecs)
#print("extended")
#print(newdecs)
lennd = len(newdecs[minHT])
deletepos = []
for y in range(1, lennd):
if y in deletepos:
continue
yp = newdecs[minHT][y]
for x in range(y + 1, lennd):
if x in deletepos:
continue
xp = newdecs[minHT][x]
if hds[yp][xp] == 0:
hds[yp][xp] = np.sum(np.logical_xor(indsv[yp].reshape(M), indsv[xp].reshape(M)))
if hds[yp][xp] < minHT:
deletepos.append(x)
print("%.2f percent" % (100.0 * y / lennd))
newdecs[minHT] = np.delete(newdecs[minHT], deletepos, axis = 0)
if len(newdecs[minHT]) <= 1:
del newdecs[minHT]
break
newdecs[minHT] = np.take(indsdec, newdecs[minHT]).tolist()
#print("deleted")
#print(newdecs)
#print(getMinimumHamming(convertIndsDecToInds(newdecs[minHT], M), M))
if len(newdecs[minHT]) == 0:
break
minHT += 2
#print("%.3f percent completed." % (i / imax * 100.0))
return newdecs
[docs]def getGoodDecsTableSmallMemory(M, K):
minHT = 4
indsiter = itertools.combinations(range(M), K)
firstivec = np.zeros(M, dtype=np.int)
firstind = np.array(next(indsiter))
firstivec[firstind] = 1
#print(firstivec)
firstdec = np.sum(np.power(2, firstind))
# Extracts the active indices having minHT >= 4
indsvec = [firstivec]
indsdec = [firstdec]
for ind in indsiter:
ivec = np.zeros(M, dtype=np.int)
npind = np.array(ind)
ivec[npind] = 1
hd = getHammingDistance(firstivec, ivec)
if hd < minHT:
continue
indsvec.append(ivec)
indsdec.append(np.sum(np.power(2, npind)))
indsvec = np.array(indsvec)
#print(np.take(indsvec, np.array([0, 1]), axis=0))
#print(len(indsvec))
#print(len(indsdec))
MCK = len(indsvec)
newdecs = {}
while True:
print("minHT = %d" % (minHT))
newdecs[minHT] = indsdec
#print(newdecs)
lennd = len(newdecs[minHT])
lstart = 0
if minHT == 4:
lstart = 1
deletepos = []
ys = np.array(list(range(lstart, lennd)))
#print(ys)
#for y in range(lstart, lennd):
yi = 0
y = ys[yi]
while True:
#if y in deletepos:
# continue
xs = np.array(list(range(y + 1, lennd)))
#print(xs)
#print(deletepos)
xs = np.setdiff1d(xs, deletepos)
if len(xs) > 0:
#print(indsvec[xs])
vxs = np.take(indsvec, xs, axis = 0)
#print(vxs.shape)
#print(vxs)
vys = np.tile(indsvec[y], len(xs)).reshape(-1, M)
#print(vys)
hds = np.sum(np.logical_xor(vxs, vys), axis = 1)
#hds = np.apply_along_axis(lambda x: getHammingDistance(indsvec[y], indsvec[x[0]]), 0, xs.reshape(1, len(xs)))
#print(hds)
#print(list(np.where(hds < minHT)[0]))
newdel = list(xs[np.where(hds < minHT)[0]])
deletepos.extend(newdel)
ys = np.setdiff1d(ys, newdel)
#print(ys)
#for x in range(y + 1, lennd):
# if x in deletepos:
# continue
# hd = np.sum(np.logical_xor(indsvec[y], indsvec[x]))
# if hd < minHT:
# deletepos.append(x)
print("%.2f percent" % (100.0 * y / lennd))
yi += 1
if yi >= len(ys):
break
y = ys[yi]
#print(deletepos)
newdecs[minHT] = list(np.delete(newdecs[minHT], deletepos, axis = 0))
if len(newdecs[minHT]) <= 1:
del newdecs[minHT]
break
if len(newdecs[minHT]) == 0:
break
minHT += 2
return newdecs
[docs]def getAllIndsBasedOnDecFile(M, K, Q):
basePath = os.path.dirname(os.path.abspath(__file__))
decfilename = basePath + "/decs/M=%d_K=%d.txt" % (M, K)
if os.path.exists(decfilename):
with open(decfilename, mode = 'r') as f:
decs = eval(f.read())
#print(decs)
print("Read " + decfilename)
minh = 0
for key in decs.keys():
if Q <= len(decs[key]):
minh = key
print(minh)
if minh > 0:
return convertIndsDecToInds(decs[minh], M)
return []
return None
[docs]def outputCPLEXModelFile(M, K, Q):
decallinds = getAllIndsBasedOnDecFile(M, K, Q)
if K > 1 and K < M-1 and decallinds == None:
print("The dec file for (%d,%d) does not exist." % (M,K))
return None
if decallinds != None and len(decallinds) > 0:
allinds = decallinds
#print(getMinimumHamming(allinds, M))
allindsvec = convertIndsToVector(allinds, M)
allindsmat = np.hstack(allindsvec).T.tolist() # MCK \times M
MCK = len(allindsmat)
else:
allinds = list(itertools.combinations(range(M), K))
allindsvec = convertIndsToVector(allinds, M)
allindsmat = np.hstack(allindsvec).T.tolist() # MCK \times M
MCK = len(allindsmat)
constraints = [" a[1] == 1;\n"]
#
#indsv = convertIndsToVector(allinds, M)
#print(getGoodDecsTable(MCK, indsdec))
#print("hds generated.")
basePath = os.path.dirname(os.path.abspath(__file__))
fname = basePath + "/inds-raw/M=%d_K=%d_Q=%d.mod" % (M, K, Q)
with open(fname, mode = 'w') as f:
f.write("int M=%d; int K=%d; int Q=%d; int MCK=%d;\n" % (M, K, Q, MCK))
f.write("int allinds[1..MCK][1..M] = " + str(allindsmat) + ";\n\n")
f.write("dvar boolean a[1..MCK];\n\n")
#
f.write("execute PARAMS {\n")
f.write(" cplex.mipemphasis = 0;\n")
f.write(" cplex.tilim = 60 * 60;\n")
f.write(" cplex.mipdisplay = 3;\n")
f.write("}\n\n")
#
f.write("minimize sum(m in 1..M) (abs(sum(q in 1..MCK)(a[q] * allinds[q][m]) - (Q * K / M)));\n\n")
#
f.write("subject to{\n")
# add constraints
f.writelines(constraints)
f.write(" sum(q in 1..MCK)(a[q]) == Q;\n")
f.write("}\n\n")
#
f.write("execute{\n")
f.write(" var f = new IloOplOutputFile(\"M=\" + M + \"_K=\"+ K + \"_Q=\" + Q + \"_obj=\" + cplex.getObjValue() + \".txt\");\n")
f.write(" f.write(a);\n")
f.write(" f.close();\n")
f.write("}\n")
print("Saved to " + fname)
return fname
[docs]def convertCPLEXOutputToInds(fname, M, K, Q):
allinds = np.array(list(itertools.combinations(range(M), K)))
decallinds = getAllIndsBasedOnDecFile(M, K, Q)
if decallinds != None and len(decallinds) > 0:
allinds = decallinds
with open(fname, mode='r') as f:
content = f.read()
content = re.sub(r'\s+', ' ', content)
content = re.sub(r'^\s+', '', content)
content = re.sub(r'\n', '', content)
content = content.replace(" ", ",")
#print(content)
inds = np.array(eval(content))
#print(inds)
#print(np.nonzero(inds)[0].tolist())
inds = np.take(allinds, np.nonzero(inds)[0], axis = 0)
return inds
[docs]def numpyToPythonStr(numpystr):
return numpystr.replace(".", "").replace(" ", " ").replace("\n\n", "\n").replace("\n ", ", ").replace(" 0", ", 0").replace(" 1", ", 1")
[docs]def convertIndsToRST(basePath, M, K, Q):
titlestr = "M = %d, K = %d, Q = %d" % (M, K, Q)
filepat = basePath + "/inds/M=%d_K=%d_Q=%d_*.txt" % (M, K, Q)
files = glob.glob(filepat)
if len(files) == 0:
print(filepat + " does not exist.")
return
fninds = files[0]
# copy the original file to the build/html/
bpath = basePath + "/../docs/build/html/db/M=%d/M=%d_K=%d_Q=%d.txt" % (M, M, K, Q)
if not os.path.exists(os.path.dirname(bpath)):
os.mkdir(os.path.dirname(bpath))
if not os.path.exists(bpath) or (os.path.exists(bpath) and os.stat(bpath).st_mtime < os.stat(fninds).st_mtime):
shutil.copyfile(fninds, bpath)
mpath = basePath + "/../docs/source/db/M=%d/" % (M)
if not os.path.exists(mpath):
os.mkdir(mpath)
fname = mpath + "M=%d_K=%d_Q=%d.rst" % (M, K, Q)
if not os.path.exists(fname) or (os.path.exists(fname) and os.stat(fname).st_mtime < os.stat(fninds).st_mtime):
with open(fname, mode = 'w') as f:
f.write("\n")
f.write("=" * len(titlestr) + "\n")
f.write(titlestr + "\n")
f.write("=" * len(titlestr) + "\n")
f.write("\n")
fn = os.path.basename(fninds)
iurl = "https://github.com/imtoolkit/imtoolkit/blob/master/imtoolkit/inds/" + fn.replace("=", "%3D")
f.write("`" + fn + " is available here. <" + iurl + ">`_\n\n")
fn = fn.replace(".txt", "")
p = Parameters(fn)
if p.Q <= 1024:
inds = np.loadtxt(fninds, dtype = np.int)
inds = inds.reshape(p.Q, p.K).tolist()
if p.Q <= 128:
at = np.array(convertIndsToMatrix(inds, M))
ats = numpyToPythonStr(str(at))
vrstr = numpyToPythonStr(str(np.array(convertIndsToVector(inds, M)).reshape(-1, p.M)))
vrstr = vrstr.replace("], ", "],\n ")
#ts = fn.replace("M=%d_"%M, "").replace("_minh=%d"%p["minh"], "").replace("_ineq=%d"%p["ineq"], "").replace("_", ", ").replace("=", " = ")
#print(ts)
#print("-" * len(ts))
f.write(".. code-block:: python\n\n")
f.write(" # minimum Hamming distance = %d\n" % p["minh"])
f.write(" # activation inequality = %d\n" % p["ineq"])
if p.Q <= 1024:
f.write(" # active indices\n")
f.write(" a = " + str(inds) + "\n")
if p.Q <= 128:
f.write(" # activation tensor\n")
f.write(" A = " + ats + "\n")
f.write(" # vector representation\n")
f.write(" " + vrstr + "\n")
else:
f.write(" # activation tensor and its vector representation are omitted.\n")
else:
f.write(" # active indices, activation tensor and its vector representation are omitted.\n")
f.write("\n")
print("The generated rst was saved to " + fname)
[docs]def main():
#print("DEBUG MODE ENABLED")
#np.set_printoptions(threshold=np.inf)
#basePath = os.path.dirname(os.path.abspath(__file__))
#convertIndsToRST(basePath, 16, 5, 2)
#quit()
if len(sys.argv) <= 1:
print("DECPARAMS_M=32")
print("DECPARAMSDO_M=32")
print("SEARCHPARAMS_M=32_P=2")
print("SEARCHPARAMSDO_M=32")
print("SEARCH_M=2_K=1_Q=2")
print("SEARCH_M=4_K=1_Q=2")
print("EVAL_dm=dic_M=2_K=1_Q=2")
print("EVAL_dm=dic_M=16_K=8_Q=16")
print("EVAL_dm=mes_M=16_K=8_Q=16")
print("EVAL_dm=wen_M=16_K=8_Q=16")
print("DECSEARCH_M=4_K=2")
print("DECEVAL_M=4_K=2")
print("DECSEARCH_M=8_K=4")
print("DECEVAL_M=8_K=4")
print("DECSEARCH_M=16_K=8")
print("DECEVAL_M=16_K=8")
print("MINH")
print("DOCMUPDATE_M=32")
print("DOCMINDEX")
print("COVERAGE")
quit()
args = sys.argv[1:]
for arg in args:
print("-" * 50)
print("arg = " + arg)
params = Parameters(arg)
basePath = os.path.dirname(os.path.abspath(__file__))
start_time = time.time()
if params.mode == "COVERAGE":
allpossibleparams = 0
hitcount = 0
M = 2
while True:
for K in range(1, M):
ps = getIMParameters(M, K)
for p in ps:
allpossibleparams += 1
M, K, Q = p[0], p[1], p[2]
fpy = glob.glob(basePath + "/inds/M=%d_K=%d_Q=%d_*.txt" % (M, K, Q))
if len(fpy) > 0:
hitcount += 1
M += 2
print("M <= %d, %d / %d = %2.2f" % (M, hitcount, allpossibleparams, 100.0 * hitcount / allpossibleparams))
if M == 32:
break
elif params.mode == "SEARCHPARAMS" or params.mode == "SEARCHPARAMSDO":
imparams = []
allpossibleparams = 0
M = 2
while True:
for K in range(1, M):
ps = getIMParameters(M, K)
for p in ps:
allpossibleparams += 1
M, K, Q = p[0], p[1], p[2]
fpy = glob.glob(basePath + "/inds/M=%d_K=%d_Q=%d_*.txt" % (M, K, Q))
if len(fpy) == 0:
imparams.append(p)
#else:
# if Q == 2 or Q * K <= M:
# print("May be wrong: /inds/M=%d_K=%d_Q=%d*.txt" % (M, K, Q))
# os.remove(fpy[0])
# fpy = glob.glob(basePath + "/decs/M=%d_K=%d.txt" % (M, K))
# if (len(fpy)) == 0:
# print("May be wrong: /inds/M=%d_K=%d_Q=%d*.txt" % (M, K, Q))
M += 2
if M > params.M:
break
print("Possible IM parameters = %d" % allpossibleparams)
print("Supported IM parameters = %d" % (allpossibleparams - len(imparams)))
print("Search coverage = %.2f percent" % (100.0 - 100.0 * len(imparams) / allpossibleparams))
imparams.sort(key = lambda x:x[2])
cmds = []
for p in imparams:
cmd = "imsearch SEARCH_M=%d_K=%d_Q=%d" % (p[0], p[1], p[2])
cmds.append(cmd)
if params.mode == "SEARCHPARAMS":
for p in range(params.P):
print("Set %d ====================================" % p)
i = 0
for cmd in cmds:
if i % params.P == p:
print(cmd)
i += 1
print("")
elif params.mode == "SEARCHPARAMSDO":
for cmd in cmds:
os.system(cmd)
elif params.mode == "DECPARAMS" or params.mode == "DECPARAMSDO":
decparams = []
allpossibleparams = 0
M = 2
while True:
for K in range(2, M - 1): # excludes K = 1 and K = M - 1
decfilename = basePath + "/decs/M=%d_K=%d.txt" % (M, K)
if not os.path.exists(decfilename):
decparams.append([M, K])
allpossibleparams += 1
M += 2
if M > params.M:
break
print("Search coverage = %.2f percent" % (100.0 - 100.0 * len(decparams) / allpossibleparams))
decparams.sort(key = lambda x:special.binom(x[0], x[1]))
cmd = "echo " + " ".join(["DECSEARCH_M=%d_K=%d" % (p[0], p[1]) for p in decparams]) + " | xargs -n1 -P15 imsearch"
if params.mode == "DECPARAMSDO":
os.system(cmd)
else:
print(cmd)
elif params.mode == "DECSEARCH":
#dectable = getGoodDecsTable(params.M, params.K)
dectable = getGoodDecsTableSmallMemory(params.M, params.K)
decfilename = basePath + "/decs/M=%d_K=%d_option=low.txt" % (params.M, params.K)
with open(decfilename, mode = 'w') as f:
f.write(str(dectable))
print("Saved to " + decfilename)
elif params.mode == "DECEVAL":
decfilename = basePath + "/decs/M=%d_K=%d.txt" % (params.M, params.K)
with open(decfilename, mode = 'r') as f:
print("Read " + decfilename)
dectable = eval(f.read())
print(dectable)
for minh in dectable.keys():
print("minh = %d" % (minh))
allinds = convertIndsDecToInds(dectable[minh], params.M)
print(allinds)
print("actual minh = %d" % (getMinimumHammingDistance(allinds, params.M)))
elif params.mode == "SEARCH":
if params.Q == 2 or params.Q * params.K <= params.M:
print("A self-evident solution is available for this setup.")
#params = Parameters("M=26_K=25_Q=16")
qstarts = np.floor(np.arange(params.Q) * (params.M - params.K) / (params.Q - 1) + 0.5)
inds = np.zeros((params.Q, params.K), dtype = np.int)
for q in range(params.Q):
inds[q] = qstarts[q] + np.arange(params.K)
outputIndsToFile(inds, params.M)
else:
fname = outputCPLEXModelFile(params.M, params.K, params.Q)
if fname != None and ".mod" in fname:
os.system("oplrun " + fname)
# Convert the obtained solution to a numpy file
fcout = glob.glob(basePath + "/inds-raw/M=%d_K=%d_Q=%d_*.txt" % (params.M, params.K, params.Q))
if len(fcout) > 0:
fcout.sort()
fname = fcout[0]
obj = 0
if "obj=" in fname:
res = re.match(r'.*_obj=(\d+)', fname)
if res:
obj = int(res.group(1))
inds = convertCPLEXOutputToInds(fname, params.M, params.K, params.Q)
outputIndsToFile(inds, params.M)
elif params.mode == "EVAL":
inds = getIndexes(params.dm, params.M, params.K, params.Q)
print(np.array(convertIndsToVector(inds, params.M)).reshape(-1, params.Q))
print("Minimum Hamming distance = %d" % getMinimumHammingDistance(inds, params.M))
print("Inequality L1 = %d" % getInequalityL1(inds, params.M))
elif params.mode == "DOCMUPDATE":
np.set_printoptions(threshold=np.inf)
M = 2
while True:
for K in range(1, M):
ps = getIMParameters(M, K)
for p in ps:
M, K, Q = p[0], p[1], p[2]
convertIndsToRST(basePath, M, K, Q)
M += 2
if M > params.M:
break
elif params.mode == "DOCMINDEX":
np.set_printoptions(threshold=np.inf)
def lfs(f):
fn = os.path.basename(f).replace(".rst", "")
p = Parameters(fn)
return (p.M * 32 + p.K) * 10000000 + p.Q
dirs = glob.glob(basePath + "/../docs/source/db/M=*")
for mdir in dirs:
print(mdir)
M = int(os.path.basename(mdir).replace("M=", ""))
titlestr = "M = %d" % M
files = glob.glob(mdir + "/M=*.rst")
files.sort(key = lfs)
fout = mdir + "/index.rst"
with open(fout, mode = 'w') as f:
f.write("\n")
f.write("=" * len(titlestr) + "\n")
f.write(titlestr + "\n")
f.write("=" * len(titlestr) + "\n")
f.write("\n\n")
f.write("This webpage provides the designed active indices in the :math:`M = %d` case.\n\n" % M)
frsts = [os.path.basename(frst).replace(".rst", "") for frst in files]
f.write(".. toctree::\n")
f.write(" :maxdepth: 2\n")
f.write(" :hidden:\n")
f.write(" \n")
for frst in frsts:
f.write(" " + frst + "\n")
f.write("\n")
for frst in frsts:
f.write("- :doc:`" + frst + "`\n")
print("The index.rst file was saved to " + fout)
elapsed_time = time.time() - start_time
print ("Elapsed time = %.10f seconds" % (elapsed_time))