# transformiert LIN(0/1). 
# Pro Spalte werden die Anzahl der positiven und negativen
# Koeffizienten gezaehlt. Dann werden in dieser Spalte Zeilen dieser Gruppe (kleiner!)
# nach oben geschoben und mit den uebrigen Zeilen darunter verrechnet, so dass in den
# unteren Spalten dann 0 drin steht. Die Anzahl der 0en pro Zeile steigt nach unten hin an.
# Dieses LIN(0/1)-Ungleichungssystem ist leichter loesbar, kann aber Loesungen zeigen,
# welche das Ursprungssystems nicht hat.
#
# J. Gamenik August 2019
import math
import random

anzZeilen = 4
anzSpalten = 4
rechtsGanzzahlig = True
random = random.Random()


# gibt Tupel von Matrix, rechter Seite und reeler Loesung zurueck
def holeBeispiel():
  global anzZeilen, anzSpalten
  #matrix = [[-1, -1, -1, 7], 
  #          [-1, 4, -1, -1], 
  #          [-1, -1, 5, -1], 
  #          [6, -1, -1, -1]
  #         ]
  #rechts = [0.1, 0.1, 0.1, 0.1]
  matrix = [[-8, 1, -7, 9], 
            [5, 7, -5, 7], 
            [-1, 10, -5, -5], 
            [5, -7, 5, -7]
           ]
  rechts = [9, -8, -1, -15]
  # Ergebnis der Reduktion
  #matrix = [[5, -7, 5, -7], [4, 3, 0, -12], [161, 0, 0, -540], [10, 0, 0, 0]]
  #rechts = [-15, -16, -884, -23]
  anzZeilen = 4
  anzSpalten = 4
  return (matrix, rechts)


# transformiere Matrix mit rechter Spalte
def transformiere(matrixO, rechtsO):
  # Neue Strukturen durch Kopie erzeugen
  matrix, rechts = kopieren(matrixO, rechtsO)
  ab = 0
  erledigt = {}
  while ab < anzZeilen:
    liste = []
    for b in range(0,anzSpalten):
      if b in erledigt:
        continue
      # ermittle Tupel (anzahl, Koeff) und gib in Liste
      anzPlus = 0
      minPlus = 1e9
      anzMinus = 0
      minMinus = -1e9
      for z in range(ab,anzZeilen):
        if matrix[z][b] > 0:
          anzPlus += 1
          if matrix[z][b] < minPlus:
            minPlus = matrix[z][b]
        if matrix[z][b] < 0:
          anzMinus += 1
          if matrix[z][b] > minMinus:
            minMinus = matrix[z][b]
      if anzPlus <= anzMinus and anzPlus > 0:
        liste.append((anzPlus,minPlus,b,minPlus))
      elif anzMinus > 0:
        liste.append((anzMinus,abs(minMinus),b,minMinus))
    if len(liste) == 0:
      break   # nur noch 0en in Spalten
    liste.sort()   # Gruppengroesse /Pivotbetrag / Spalte
    #print(liste)
    # nehme die Spalte mit einer moeglichst kleinen Gruppe und kleinem
    # Koeffizienten,. ggf werden Zeilen getauscht
    plus = True if liste[0][3] > 0 else False
    spalte = liste[0][2]
    optKoeff = liste[0][1]
    # alle Zeilen dieser Gruppe nach oben bringen
    for lauf in range(0,liste[0][0]):
      if ((matrix[ab+lauf][spalte] <= 0 and plus) or (matrix[ab+lauf][spalte] >= 0 and not plus)):
        for lauf2 in range(ab+lauf+1,anzZeilen):
          if ((matrix[lauf2][spalte] > 0 and plus) or (matrix[lauf2][spalte] < 0 and not plus)):
            vertauscheZeile(ab+lauf, lauf2,matrix,rechts)
            break
    # Zeile mit optKoeff sollte oben stehen [ab]
    if abs(matrix[ab][spalte]) != optKoeff:
      for lauf in range(ab+1,anzZeilen):
        if abs(matrix[lauf][spalte]) == optKoeff:
          vertauscheZeile(ab, lauf,matrix,rechts)
          break 
    # alle restlichen Zeilen der anderen Gruppe mit oberstem verrechnen
    abV = ab + liste[0][0]
    for lauf in range(abV,anzZeilen):
      if ab != lauf and matrix[lauf][spalte] != 0:
        verrechneZeile(ab,lauf,spalte,matrix,rechts)
    ab = ab + liste[0][0]
    erledigt[spalte] = 1
  return (matrix, rechts, erledigt)


def kopieren(matrixO, rechtsO):
  matrix = []
  for a in range(0,anzZeilen):
    zeile = []
    for b in range(0,anzSpalten):
      zeile.append(matrixO[a][b])
    matrix.append(zeile)
  rechts = []
  for a in range(0,anzZeilen):
    rechts.append(rechtsO[a])
  return (matrix, rechts)


def vertauscheZeile(a,zeile,matrix,rechts):
  if (a != zeile):
    for b in range(0,anzSpalten):
      sicher = matrix[a][b]
      matrix[a][b] = matrix[zeile][b]
      matrix[zeile][b] = sicher
    sicher = rechts[a]
    rechts[a] = rechts[zeile]
    rechts[zeile] = sicher


def verrechneZeile(z1,z2,spalte,matrix,rechts):
  pivot = matrix[z1][spalte]
  vergleich = matrix[z2][spalte]
  durch = ggT(pivot, vergleich)
  #print("Pivot {0} und Vergleich {1}".format(pivot,vergleich))
  #print(matrix[z1])
  #print(matrix[z2])
  if pivot * vergleich >= 0:
    raise RuntimeError("Verrechnung nicht moeglich")
  pivot = abs(pivot)
  vergleich = abs(vergleich) 
  for a in range(0,anzSpalten):
    matrix[z2][a] = matrix[z1][a] * vergleich + matrix[z2][a] * pivot
    matrix[z2][a] //= durch
  if matrix[z2][spalte] != 0:
    print(matrix[z2])
    raise RuntimeError("Verrechnung hat nicht geklappt")
  rechts[z2] = rechts[z1] * vergleich + rechts[z2] * pivot
  if rechts[z2] % durch == 0:
    rechts[z2] //= durch
  else:
    rechts[z2] /= durch
  #print("Danach", end="")
  #print(matrix[z2], end="")
  #print(rechts[z2])
  # Optimierung: wenn alle Koeffizienten einer Zeile ggT > 1 haben, dann teilen
  ggt = matrix[z2][0]
  for a in range(1,anzSpalten):
    ggt = ggT(matrix[z2][a],ggt)
  if ggt > 1:
    #print("Alle Koeffizienten haben ggT " + str(ggt))
    ggt2 = ggT(ggt,rechts[z2])
    if (rechtsGanzzahlig and ggt2 > 1) or rechtsGanzzahlig == False:
      if rechtsGanzzahlig:
        ggt = ggt2
      #print("Teilen durch ggT " + str(ggt))
      for a in range(0,anzSpalten):
        matrix[z2][a] //= ggt
      if rechts[z2] % ggt == 0:
        rechts[z2] //= ggt
      else:
        rechts[z2] /= ggt


def ggT(w1, w2):
  w1 = abs(w1)
  w2 = abs(w2)
  if w1 < w2:
    hilfe = w1
    w1 = w2
    w2 = hilfe
  while (w2 > 0):
    rest = w1 % w2
    w1 = w2
    w2 = rest
  return w1


if __name__ == "__main__":
  matrix, rechts = holeBeispiel()
  matrix2, rechts2 = transformiere(matrix, rechts)
  print(matrix)
  print (rechts)