gap for projections of graphic matroids

projection.sage

@@ -0,0 +1,55 @@
+# cogirth-packing gap of projections of graphic matroids
+
+from sage.all import *
+from sage.matroids.all import *
+from sage.graphs.all import *
+import gurobipy as gp
+from gurobipy import GRB
+
+env = gp.Env(empty=True)
+env.setParam("OutputFlag",0)
+env.start()
+
+
+def cogirthip(bases, integral=true):
+    model = gp.Model("mip1",env=env)
+    # model.Params.LogToConsole = 0
+    groundset=frozenset()
+    for B in bases: groundset=groundset|frozenset(B)
+    x = dict()
+    if integral:
+        for e in groundset: x[e]=model.addVar(vtype=GRB.BINARY)
+    else:
+        for e in groundset: x[e]=model.addVar(vtype=GRB.CONTINUOUS,lb=0)
+    model.setObjective(gp.quicksum([x[e] for e in groundset]), GRB.MINIMIZE)
+    for B in bases:
+        model.addConstr(gp.quicksum([x[e] for e in B])>=1)
+    model.optimize()
+    return model.ObjVal
+
+
+cnt=0   # actual number of instances tested
+f = lambda g: g.is_connected()
+for n in range(4,5):
+    for g in filter(f, graphs(n)):
+        g=graphs.CompleteGraph(4)
+        A=g.incidence_matrix()
+        n,m = A.dimensions()
+        # enumerate all vectors in F_2^n with even number of 1s
+        V = VectorSpace(GF(2), n)
+        for v in filter(lambda v: v.hamming_weight() % 2 == 0 and v!=0, V):
+            v_col=matrix(v).transpose()
+            A_t=A.augment(v_col)
+            # print(A_t)
+            M=Matroid(matrix=A_t,field=GF(2))/m     #contract the last element
+            # print(M)
+            bases=M.bases()
+            strength=cogirthip(bases,integral=false)
+            cogirth =cogirthip(bases,integral=true)
+            gap = cogirth/strength
+            if gap > 2.01:
+                print(f"bad example! gap={gap}")
+                with open("projection.out", "a") as f:
+                    f.write(str(A_t)+"\n\n")
+            else:
+                print(f"n={n},m={m},gap={gap}")