// GNFS - General Number Field Sieve
// Versucht, eine einfache Programmierung des verallg. Zahlkörpersiebs zu machen.
// Siehe auch Beschreibung von "Matthew E. Briggs", Kapitel 5
// Beginn: Sa 17.03.12
// Autor: J. Gamenik


// Anm.: Parameter sind Zahl und Basis (g-adische Zahl), im folgenden n und m bezeichnet
// Z[alpha] ist der algebraische Zahlkörper, d.h. m ist ganzzahlige Nullstelle des Polynom mod n,
// alpha eine Nullstelle aus dem komplexen Bereich C des gewählten Polynoms.
// Das Polynom soll über Z irreduzibel (sonst sehr einfach faktorisierbar) sowie das Minimalpolynom zu alpha sein,
// d.h. es gibt kein Polynom kleineres Grades mit alpha als Nullstelle, d.h. alpha kommt nicht mehrfach vor.
// Die Norm und die Abbildung (Ringhomomorphismus nach Z) von (a + b * alpha) ist errechenbar, nämlich
// (-b)^d * f(-a/b).
// Es werden nun simultan für (a,b) die Produkte (a + b * m) sowie N(a + b * alpha) als Quadrate dargestellt,
// wobei modulo n gerechnet wird.
// Im zweiten Ausdruck ist dann aber nicht sicher, dass die Produkte (a + b * alpha) dann auch wirklich ein
// Quadrat in Z[alpha] ergeben, wenn das Produkt der Normen ein Quadrat ist !

// Anm.: Im Spezialfall Grad d=2 wird hier eine Vollsuche nach dem Wurzelpolynom gemacht, was nicht
// praxistauglich ist.

import java.math.BigInteger;
import java.util.BitSet;
import java.util.List;
import java.util.ArrayList;
import java.io.IOException;

public class GNFS1 extends Thread
{
  private static long n;   // zu faktorisierende Zahl
  private static long[] poly;  // zu ben. Polynom
  private long m;  // Zahlenbasis
  private GNFS_Paar [] liste;  // Liste mit Relationen
  private static int d;  // Grad des Polynoms
  private final long E;  // Grenze für a und b
  private final long B1;  // B1-glatt für Zahl
  private final long B2;  // B2-glatt für Norm vom Polynom
  public final static long primBasis [] = {-1,2,3,5,7,11,13,17,19,23,29}; // ,31,37,41,43};
  private long x, beta;  // Ergebnis des Verrechnungsschrittes

  private boolean nurPlus = true;  // Einfluss auf Polynomgenerierung
  private final static boolean DEBUG = false;
  private final static boolean DEBUG2 = false;  // zerlegungen
  public final static boolean DEBUG3 = false;  // Verrechnungen; Achtung: in Hilfsklasse gibt es noch ein DEBUG3
  public final static boolean DEBUG4 = false;  // nicht geglückte Zerlegungen

  public static void main (String [] args)
  {
    long eing1 = 45113, eing2 = 31;  // Briggs-Beispiel n=45113 mit m=31
    if (args.length >= 2) eing2 = Integer.parseInt(args[1]);
    if (args.length >= 1) eing1 = Integer.parseInt(args[0]);
    System.out.println ("Achtung: Experimenteller Status, noch nicht zu Ende programmiert!");
    if (eing1 < 1 || eing2 < 1 || eing1 < eing2)
    {
      System.out.println ("Syntax: java GNFS1 <zahl> <basis> zahl=" + eing1 + ", basis=" + eing2);
      System.exit(1);
    }

    GNFS1 sieb = new GNFS1(eing1, eing2);
    sieb.start();
  }

  public GNFS1 (long p_zahl, long p_basis)
  {
    if ((p_zahl < 1) || (p_basis < 2)) { System.out.println ("Ungültige Argumente"); System.exit(1); }
    n = p_zahl; m = p_basis; 
    E = 37; B1 = 29; B2 = 29;  // Schranken, um über Zerlegbarkeit entscheiden zu können
    if (primBasis[primBasis.length-1] != B1) { System.out.println ("Logikfehler: Primbasis veraltet"); System.exit(1); }
    liste = new GNFS_Paar[primBasis.length * 2];  // da zwei Ausdrücke korresp. zu zerlegen sind
  }

  public static long holeN() { return n; }  // n private und static, soll nicht verändert werden aussen
  public static long [] holePoly() { return poly; }
  public static int holeGrad() { return d; }

  // Hauptschleife
  public void run()
  {
    System.out.println ("--> Zahl: " + n + ", g-adische-Basis " + m);
    if (DEBUG) System.out.println ("Schranken: E=" + E + ", B1=" + B1 + ", B2=" + B2); 
    System.out.println ("Notation: Ergebnis [(a + b * m) , N(a + b * alpha)]\n---------------------------------------------------");

    // 1. Schritt: Polynomwahl
    polynomwahl_machen();
    // 2. Schritt: Sieb
    // 3. Schritt: Gleichungssystem aufbauen
    // 4. Schritt: Wurzelberechnung, wird innerhalb des Gls-Aufbau gemacht
    sieben_aufbauen();
  }


  // -------------------------------------------------
  // hier die Ausprogrammierung der einzelnen Schritte
  // -------------------------------------------------

  // I) bestimmt das Polynom anhand der g-adischen Darstellung der Zahl n zur Basis Var. m
  private void polynomwahl_machen()
  {
    // wähle monisches, irreduzibles Polynom vom Grad d zur Basis m, wo m eine Nullstelle modulo n bildet
    // m wurde schon im Konstruktor festgelegt.  m^d ~ n, d = (int) (ln(n) / ln(m))
    d = (int) ( Math.log(n) / Math.log(m) );
    if (d <= 1)
    {
      System.out.println ("Polynomgrad <= 1, nicht geeignet, Basis kleiner wählen");
      System.exit(1);
    }
    long feld [] = null;
    // Koeffizienten sind nur >= 0, also genau die m-äre Darstellung
    feld = new long[d+2 /* 2 weil ggf. um 1 erhöht */]; long zerl = n;
    for (int i = 0; i <= d; ++i)
    {
      feld[i] = zerl % m; zerl = zerl / m;
    }
    if (feld[0] == 0) System.out.println ("... Hinweis: ganzzahliger Teiler " + m);
    // Normieren
    if (feld[d] != 1)
    {
      feld[d-1] += (feld[d]-1) * m; feld[d] = 1;
      // Anm.: Normierung scheint wichtig zu sein, weil sonst es bei der Suche nach Wurzel-
      // Polynomen auch bei syst. Komplettsuche keinen Treffer geben muss. Grund ? z.B.
      // n=19, m=3 oder n=53, m=5
    }

    if (!nurPlus)
    {
      // Koeffizienten können auch < 0 sein, z.B. 15 = x^4 - 1 bei m=2
      long mitte = m / 2; boolean anz = false;
      for (int i = d; i >= 0; --i)
      {
        if (feld[i] > mitte)
        {
          // m = 3: 25 = 2|2|1  --> 1|0|0|-2
          for (int j = 0; j <= i; ++j) { feld[j] = feld[j] - m; feld[j+1]++; } anz = true;
        }
        else if (feld[i] < -mitte)
        {
          // m = 3: -8 = -2|-2  --> -1|0|1
          for (int j = 0; j <= i; ++j) { feld[j] = m + feld[j]; feld[j+1]--; } anz = true;
        }
        if (anz)
        { 
          if ((d+1) < feld.length && feld[d+1] != 0) ++d;
          if (DEBUG) { for (int k = d; k >= 0; --k) { System.out.print(feld[k]+","); } System.out.println(""); }
          anz = false; i = d+1 /* wieder von oben los */;
        }
      }
      // zu tun: Normierung sicherstellen
    }

    // Test, ob Polynom passt ...
    long prf = feld[d];
    for (int i = d-1; i >= 0; --i) { prf = prf * m + feld[i]; }
    if (prf != n) { System.out.println ("Unterschied: " + prf + " und " + n); System.exit(1); }

    // Test auf Irreduzibilitaet entfällt momentan (sonst: ist faktorisierbar)
    // Test auf tatsächlich Minimalpolynom von alpha (Nullstelle von Polynom im komplexen) entfällt,
    // d.h. die Nullstellen müssen nicht einfach sein, sondern, können mehrfach vorkommen.

    poly = feld;
    ausgabe_polynom();
  }  

  // Gibt das Polynom zur Basis m aus
  private void ausgabe_polynom()
  {
    System.out.println("Polynom: ");
    boolean ausg = false;
    for (int k = d; k >= 0; --k) 
    {
      if(poly[k] != 0) 
      {
        if (ausg) System.out.print(" + ");
        System.out.print(poly[k]); 
        if (k > 0) System.out.print("*" + m + "^" + k); ausg=true;
      }
    }
    System.out.println("");
  }


  // II) führt den Siebschritt aus. Sucht Paare a und b, so dass zwei davon abgel. Ausdrücke glatt sind.
  // III) baut da Gleichungssystem auf
  private void sieben_aufbauen()
  {
    int erfolg = 0, misserfolg = 0;
    // Anm.: In Z betrachtet, also nicht mod. n gerechnet, scheint der Funktionswert der Norm mit steigenden b und innerhalb
    // dessen mit steigenden (absolut) a anzuwachsen. Es treten oft Primzahlen auf (Primideale).
    for (long b = 1; b <= E; ++b)
    {
      // Beachte Symmetrie a, b und -a,-b ist gleich bzl. Zerlegbarkeit in Bereich Z, also nur a ins neg. laufen lassen
      // Muster: 1, -1, 2, -2, 3, -3 u.s.w.
      long a = +1; boolean wechsel = true;
      // for (long a = -E; a <= E; ++a)
      do
      {
        boolean testen = true;
        if ((a == 0) || (a % m == 0)) testen = false;  // bringt nichts
        // ggT(a,b) soll 1 sein
        if (! BigInteger.valueOf(a).gcd(BigInteger.valueOf(b)).equals(BigInteger.ONE)) testen = false;
        if (!testen)
        {
          // a noch setzen, siehe auch Schleifenende
          if (wechsel) { a = -a; wechsel = false; } // 1 --> -1
          else if (a < 0) { --a; wechsel = true; } // -1 --> -2
          else if (a > 0) { ++a; wechsel = true; } // 1 --> 2
          continue;
        }

        // Bed. 1: a + b*m ist B1-glatt
        long ausd1 = (a + b * m) % n;
        GNFS_Teil_Ganzzahl teil1 = zerlegen (ausd1, B1);

        // Bed. 2: N(a + b*alpha) = (-b)^d * f(-a/b) ist B2-glatt
        // Um hier die Norm möglichst oft teilen zu können, ist offenbar folgendes zu tun: a/b verringern, m vergrößern
        long ausd2 = norm(a, b);
        GNFS_Teil_Ganzzahl teil2 = zerlegen(ausd2, B2);

        if (ausd1 == 0 || ausd2 == 0)
        {
          // Beispiel: n=11, m=3, x^2 + 2, a=5, b=2 erzeugt 0 und 0. Aber 2x + 5 kein Quadrat % (x^2 + 2)
          System.out.println ("  Sonderfall: a=" + a + " und b=" + b + " --> [" + ausd1 + " , " + ausd2 + "]");
        }
        else if (teil1 != null && teil2 != null)
        {
          System.out.println ("  Möglich: a=" + a + " und b=" + b + " --> [" + ausd1 + " , " + ausd2 + "]");
          if (DEBUG2) { teil1.ausgabe(); teil2.ausgabe(); }
          boolean fertig = einfuegen(teil1, teil2, a, b);
          if (fertig) { System.out.println ("... Fertig"); return; }
        }
        else if (DEBUG4)
        {
          System.out.println ("  Nicht zerlegbar: a=" + a + " und b=" + b + " --> [" + ausd1 + " , " + ausd2 + "]");
        }

        if ((teil1 != null) && (teil2 != null)) ++erfolg;
        else ++misserfolg;

        // a aktualisieren
        if (wechsel) { a = -a; wechsel = false; } // 1 --> -1
        else if (a < 0) { --a; wechsel = true; } // -1 --> -2
        else if (a > 0) { ++a; wechsel = true; } // 1 --> 2
      } 
      while ((a > -E) && (a < E)); // Ende a
    } // Ende b

    System.out.println ("Statistik: Erfolg=" + erfolg + ", Misserfolg=" + misserfolg);
  }


  // zerlegt den Ausdruck numerisch über der Faktorbasis. Hier ist eine Grenze zu beachten
  // Testet ausd und -ausd mod. modulo
  private GNFS_Teil_Ganzzahl zerlegen (long ausd, long grenze)
  {
    GNFS_Teil_Ganzzahl erg = new GNFS_Teil_Ganzzahl();
    if (ausd == 0) return erg; // Endlosschleife vermeiden
    long sicher = ausd;

    if (ausd == 31 && n == 41)
    {
      int stopp = 1;  // f. Debugger
    }

    for (short runde = 0; runde < 2; ++runde)  // 0.te Runde ausd, 1.te Runde modulo - ausd
    {
      if (ausd < 0) { ausd = -ausd; erg.hoch.set(0); } // -1 an primBasis[0]

      for (int i = 1; i < primBasis.length && ausd != 1; ++i)
      {
        if (primBasis[i] > grenze) break;
        long anz = 0;
        while ((ausd % primBasis[i]) == 0) 
        { 
          ausd /= primBasis[i]; anz++;
          if (anz == 2) { erg.wurz = (erg.wurz * primBasis[i]) % n; anz = 0; }
        }
        if (anz == 1) erg.hoch.set(i);
      }
      // ist der Rest ein Quadrat, kann man noch was machen
      if (ausd != 1)
      {
        long wrz = ist_wurzel (ausd);   // (long) Math.sqrt(ausd);
        if (wrz >= 0)
        { 
          // System.out.println ("Glücksfall: Rest " + ausd + " ist Quadrat"); 
          erg.wurz = (erg.wurz * wrz) % n;
          return erg;  // Fehlt noch Übergabe der Wurzel aus Rest, nicht in Faktor-Basis
        }
      }
      // if (ausd != 1) System.out.println ("Rest " + ausd);

      if (ausd == 1) return erg; // sind alle Faktoren raus ?
      if (runde == 0) 
      { 
        // abhängig, ob Ausgangswert pos. oder neg. (9 --> -(n-9), -9 --> (n+9) )
        if (sicher > 0)
          ausd = -(n - sicher); // minus vor Klammer hier, damit oben [0] gesetzt
        else
          ausd = (n + sicher);
        erg = new GNFS_Teil_Ganzzahl();  // wiederherstellen
      }
    }

    return null;
  }


  // Berechnet die Norm bzgl. den Polynomausdrucks (a + b * alpha), siehe oben; a, b > 0
  // Führt das ganze symbolisch als Bruch durch. 
  private long norm (long a, long b)
  {
    long zaehler = poly[d]; long nenner = 1;
    for (int i = d-1; i >= 0; --i)
    {
      long nNeu = b*nenner; zaehler = (zaehler * (-a)) + (poly[i] * nNeu); nenner = nNeu;
      // System.out.println ("Zaehler " + zaehler + ", Nenner " + nenner);
    }
    // Idee: kann man hier modulo n rechnen ? Ich glaube ja, weil dann zu 11303 2 mit a=-2, b=1
    // der Funktionswert a + b * m und N(a + b * alpha) = 0 ist
    zaehler = zaehler % n;  // ggf. schon oben verringern
    return (((d & 1) == 1) ? -zaehler: zaehler );
  }

  // Fügt die Teilergebnisse in die Matrix ein, und macht Verrechnungen. 
  // 3 mögliche Folgen:
  // a) falsch, wenn triviales Paar erzeugt (Faktor 1, keine Exponenten)
  // b) wahr, wenn nicht triviales Paar erzeugt (Wurzel, aber keine Exponenten)
  // c) hinzugefügt, ggf. mit Verrechnungen
  private boolean einfuegen (GNFS_Teil_Ganzzahl teil1, GNFS_Teil_Ganzzahl teil2, long a, long b)
  {
    GNFS_Paar neu = new GNFS_Paar(teil1, teil2, a, b);

    int i; int hNeu = neu.hoechstens(); if (DEBUG) System.out.println ("Höchster Index neu ... " + hNeu);

    // if (liste[0] == null) { liste[0] = neu; neu.ausgabe(); return (hNeu < 0); }  // Sond. noch leer
    // Achtung: Verzichte erstmal auf Ordnung der Liste
    for (i = 0; i < liste.length; ++i)
    {
      if (liste[i] == null) break;  // alles durchsucht
      int hBisher = liste[i].hoechstens();
      if (hBisher == hNeu)
      {
        // damit verrechnen
        if (DEBUG3) System.out.println ("Verrechne neu mit Gleichung an " + i);
        GNFS_Paar erzeugt = liste[i].verrechnen(neu);
        int hNeu2 = erzeugt.hoechstens(); if (DEBUG) System.out.println ("Höchster Index nach verrechnet ... " + hNeu2);
        if (hNeu2 == hNeu)
        {
          System.out.println ("Fehler: Grad bisher " + hNeu + ", jetzt " + hNeu2);
        }
        if ((hNeu2 < 0) && (erzeugt.wurzel_liefern(false) == 1) && (  erzeugt.wurzel_liefern(true) == 1)) 
        {
          System.out.println ("Abbruch, weil Identität erzeugt");  //M Fall a)
          return false;  // trivial, liefert keine Info mehr
        }
        if (hNeu2 < 0)
        {
          // über Faktorbasis zerlegt
          System.out.print ("--> Gefunden: Polynom Z[alpha]/f(x)%n, ");
          hNeu = hNeu2; neu = erzeugt; break;
        }
        hNeu = hNeu2; neu = erzeugt; i = -1; // dann mit dem erzeugten wieder von vorne los
      }
    }


    if (hNeu < 0)
    {
      // die beiden Ausdrücke sind komplett glatt. Faktorisierung versuchen ...
      /* (S. 18)
        Unglücklicherweise liegt hier nur eine notwendige, aber keine hinreichende Bedingung
        dafür vor, dass das Produkt (a  +  b * alpha) ein Quadrat in Z[alpha] ist, wenn das
        Produkt der Normen ein Quadrat ist. Der Grund ist, dass Normabbildung nicht injektiv ist,
        d.h. es ist möglich, ein Produkt mit quadr. Norm zu konstruieren, ohne dass dieses
        Produkt selbst ein Quadrat in Z[alpha]. --> Hier trennen sich SNFS und GNFS.
        SNFS setzt einen Hauptidealring voraus.
      */

      // Ergebniswerte übernehmen (Normen).
      this.x = neu.wurzel_liefern(false); this.beta = neu.wurzel_liefern(true);
      System.out.println ("x ist " + this.x + ", Norm(beta) ist " + this.beta);

      // aus den gespeicherten Paaren (a,b) das (quadratische) Polynom berechnen
      long [] erg = neu.polynom_liefern();
      System.out.print("Quadrat-Polynom ermittelt: "); boolean ausg = false;
      for (int j = erg.length-1; j >= 0; --j) 
      {
        if (erg[j] != 0) { 
          if (ausg) System.out.print(" + "); 
          System.out.print(erg[j]); 
          if (j > 0) System.out.print("*x^" + j + " "); ausg = true; }
       }
       System.out.println("");

       // jetzt wirds schwierig: aus diesem Polynom eine Wurzel in Z[alpha] finden, notiert beta.
       long [] wurz = wurzel_polynom_liefern (erg);

       long y = -1;
       if (wurz != null)
       {
         // Ganzzahl y bestimmen, indem m in das Wurzelpolynom eingesetzt wird.
         y = abbildung_liefern(wurz);
         System.out.println ("Abgebildeter Wert y ist " + y);
       }
       else
       {
         y = this.beta;  // als Ersatz
       }

       // ggT aus (x +- y) betrachten
       for (short runde2 = 0; runde2 < 2; ++runde2)
       {
         long diff1 = ((runde2 == 0) ? (this.x - y) : (this.x + y)) % n;
         if (diff1 == 0) System.out.println ("Sind gleich");
         else {
           BigInteger ggt1 = BigInteger.valueOf(diff1).gcd(BigInteger.valueOf(n));
           if(! ggt1.equals(BigInteger.ONE) )
           {
             System.out.println ("Treffer " + ggt1); return true;
           }
         }
       }

       System.out.println(""); warten();
       return false /* true */; // Fall b)
    }

    if ((i < 0) || (i >= liste.length)) 
    {
      System.out.println ("Fehler: Einfügen an " + i + " von " + liste.length); 
      System.exit(1);
    }

    // Sonst Fall c): füge neu an der ersten freien Position ein (nicht sortiert !)
    if ((hNeu >= 0) || (neu.wurzel_liefern(false) != 1) || (neu.wurzel_liefern(true) != 1))
    {
      if (DEBUG) System.out.println ("Neu an " + i + " eingefügt:");
      liste[i] = neu;
      neu.ausgabe();
    }

    return false;
  }

  // wartet auf Benutzereingabe
  public static void warten()
  {
    // if (true) return;
    try
    {
      int zeichen;
      while ((zeichen = System.in.read()) != '\n');
    } catch (IOException ex) {}
  }

  // Wurzel modulo möglichst berechnen, also ggf. Vielfache von modulo addieren
  long ist_wurzel (long ausd)
  {
    while (ausd < 0) ausd += n;  // danach: >= 0
    long wert = ausd;
 
    if (wert <= 1) return wert;
    long weg = 1; 
    // 4er raus
    while ((wert & 1) == 0)
    {
      if ((wert & 3) == 0) { weg <<= 1; wert >>= 2; }
      else if ((wert & 1) == 0) wert += n;
    }
/*
    long [] quWeg = {9, 25, 49, 121, 169, 289};  // bis 17^2
    long [] wWeg = {3, 5, 7, 11, 13, 17};
    do
    {
      if ((wert & 3) == 3) wert += (2 * n);
      if ((wert & 7) == 5) wert += (4 * n);
      // ist jetzt 1 mod 8, Rauswurf ungerade Quadrate ändert daran nichts, also
      boolean reduziert = false;
      for (int i = 0; i < quWeg.length; ++i)
      {
        if ((wert % quWeg[i]) == 0) { wert /= quWeg[i]; weg = (weg * wWeg[i]) % n; 
              System.out.println("Weg " + weg); reduziert = true;}
      }
      if (!reduziert) break;
    } 
    while (true);
*/

    long test = (long) Math.sqrt(wert);
    if ((test * test) == wert) 
    {
      // System.out.println ("... Wurzel gefunden " + weg + ", " + test + "," + ausd);
      long teil = ((weg * test) % n);
      if ( ((teil * teil) - ausd) % n != 0)
      {
        System.out.println ("Wurzelberechnung ist falsch"); System.exit(1);
      }
      return ((weg * test) % n);
    }

    return -1;  // nichts gefunden
  }


  // versucht, eine Wurzel in Z[alpha] zu finden (schwer, arbeitet hier händisch)
  private long [] wurzel_polynom_liefern(long [] quadrat)
  {
    if (quadrat == null || quadrat.length == 0) { System.out.println ("Leeres Quadratpolynom"); System.exit(1); }
    boolean komplett = true;

    // Schaue Spezialfälle an, z.B. Quadratpolynom vom Grad 0
    int h_grad = quadrat.length-1;
    for (; h_grad >= 0; --h_grad) if (quadrat[h_grad] != 0) break;

    if (!komplett)
    {

    // Sonderfall Polynom vom Grad 2: dann hier versuchen zu knacken
    if ((h_grad == 1) && (d == 2))
    {
      for (int i = 0; i <= 2; i++) quadrat[i] = (quadrat[i] + poly[i]) % n;
      h_grad = 2;
    }
    if (h_grad == 2); //System.out.println ("Sonderfall Polynom Grad 2");
    else if (h_grad != 0) return null;  // sonstige ungeraden Polynome

    // Versuche, Wurzeln zu ziehen
    if (h_grad == 2)
    {
      long koeff2 = ist_wurzel(quadrat[2]); if (koeff2 < 0) return null;
      long koeff1 = ist_wurzel(quadrat[0]); if (koeff1 < 0) return null;
      // noch mittleres Element testen
      long vergleich = (2 * koeff1 * koeff2) % n;
      if ((vergleich - quadrat[1]) % n == 0)
      {
        long [] rck = new long[2]; rck[0] = koeff1; rck[1] = koeff2; return rck;
      }
      else return null;
    }
    else
    {
      long teilerg = ist_wurzel(quadrat[0]);
      if (teilerg != -1)
      {
        long [] rck = new long[1]; rck[0] = teilerg; return rck;
      }
    }

    } // Ende nicht komplett
    else if (d == 2)
    {  // Start komplett

      // Ann.: poly[2] != 0
      long durch = 0;  // (1 / höchste Polynomkoeff.)
      try
      {
         durch = (poly[2] == 1) ? 1 : BigInteger.valueOf(poly[2]).modInverse(BigInteger.valueOf(n)).longValue();
      } 
      catch (ArithmeticException ex) { System.out.println ("ggT gefunden"); System.exit(0); }

      System.out.println ("Vollsuche, weil Basispolynom Grad 2 hat; kann dauern ...");
      // alle Möglichkeiten durchtesten, nicht praxistauglich bei gr. Zahlen
      // Teste alle Polynome der Form (lauf1*x + lauf2)^2
      long bisher = System.currentTimeMillis();
      for (long lauf1 = 0; lauf1 < n; ++lauf1)  // Grad 1
      {
        for (long lauf2 = 0; lauf2 < n; ++lauf2) // Grad 0
        {
          long k2 = (lauf1 * lauf1) % n; 
          long k1 = (2 * lauf1 * lauf2) % n; 
          long k0 = (lauf2 * lauf2) % n;
          if (k2 != 0)
          {
            // Ganzzahliges abziehen
            long weg = (k2 * durch) % n;
            k2 = (k2 - weg * poly[2]) % n;
            k1 = (k1 - weg * poly[1]) % n;
            k0 = (k0 - weg * poly[0]) % n;
            if (k2 != 0) { System.out.println ("Fehler beim Polynomerzeugung"); System.exit(1); }
          }

          // vergleichen, alle Koeffizienten können auch neg. sein, bei Vz-Ungleichheit handeln.
          if (k2 < 0 && quadrat[2] > 0) k2 += n;
          else if (k2 > 0 && quadrat[2] < 0) k2 -= n;
          if (k1 < 0 && quadrat[1] > 0) k1 += n;
          else if (k1 > 0 && quadrat[1] < 0) k1 -= n;
          if (k0 < 0 && quadrat[0] > 0) k0 += n;
          else if (k0 > 0 && quadrat[0] < 0) k0 -= n;
          // System.out.println ("Polynom " + k2 + " , " + k1 + " , " + k0);

          if ((k2 == quadrat[2]) && (k1 == quadrat[1]) && (k0 == quadrat[0]))
          {
            System.out.print ("Treffer Komplettsuche, Polynom ");
            long rueck [] = new long [2];
            rueck[1] = lauf1; rueck[0] = lauf2; 
            System.out.println(lauf1 + "*x + " + lauf2);
            return rueck;
          }
        }
        if ((lauf1 % 1000) == 999)
        {
          long diff = System.currentTimeMillis() - bisher;
          long gesamt = (diff / lauf1) * n;
          long rest = gesamt - diff;
          System.out.println("Noch max sek: " + (long) (rest / 1000.0) );
        }
      }
      System.out.println ("Ende Vollsuche");      
    } // Ende komplett

    return null;
  }

  // Bildet einen Wert auf Z mod n ab. Horner-Schema
  private long abbildung_liefern (long [] eingabe)
  {
    if (eingabe == null || eingabe.length == 0) { System.out.println ("Leeres Wurzelpolynom"); System.exit(1); }

    long rueck = eingabe[eingabe.length-1];
    for (int i = eingabe.length-2; i >= 0; --i)
    {
      rueck = (((rueck * m) % n) + eingabe[i]) % n;
    }

    return rueck;
  }

}



// Hilfsklasse: 
// speichert zu einem Produkt die Wurzeln und binären Exponenten, d.h. alles aus den ganzen Zahlen, bei (a + b * alpha) ist das die Norm
class GNFS_Teil_Ganzzahl
{
  public long wurz;  // der Wurzel-Teil des quadr. Faktors
  public BitSet hoch;  // Stelle, wo ein Faktor nicht quadratisch, sondern nur einfach vorkommt
                       // Beispiel: Prod. 3 * 5^2 * 11 hat wurz = 5 und Index 2 f. prim 3 und Ind. 5 für prim 11

  // Standard-Konstruktor
  public GNFS_Teil_Ganzzahl ()
  {
    wurz = 1; hoch = new BitSet(GNFS1.primBasis.length);
  }
  // Ausgabe:
  public void ausgabe()
  {
    if (wurz != 1) System.out.print ("Wurzel " + wurz + ", ");
    if (hoch.isEmpty() && wurz != 1) System.out.println("");
    else if (!hoch.isEmpty())
    {
      System.out.print("einfache Fakt. "); 
      for (int i = 0; i < hoch.length(); ++i) if (hoch.get(i)) System.out.print(GNFS1.primBasis[i]+",");
      System.out.println(""); 
    }
  }

  // Verrechnen
  public GNFS_Teil_Ganzzahl verrechnen (GNFS_Teil_Ganzzahl mit)
  {
    if (GNFS1.DEBUG3) System.out.println (wurz + " und " + hoch);
    if (GNFS1.DEBUG3) System.out.println (mit.wurz + " und " + mit.hoch);
    GNFS_Teil_Ganzzahl erg = new GNFS_Teil_Ganzzahl();
    erg.wurz = (this.wurz * mit.wurz) % GNFS1.holeN();
    if (this.hoch.get(0) != mit.hoch.get(0)) erg.hoch.set(0);  // Vz
    for (int i = 1; i < GNFS1.primBasis.length; ++i)
    {
      if (this.hoch.get(i) && mit.hoch.get(i)) { erg.wurz = (erg.wurz * GNFS1.primBasis[i]) % GNFS1.holeN(); } 
      else if (this.hoch.get(i) != mit.hoch.get(i)) erg.hoch.set(i);
      // im Fall beide gesetzt oder beide leer ist neuer Exponent 0, das ist Ziel !
    }
    if (GNFS1.DEBUG3) System.out.println ("... ergibt " + erg.wurz + " und " + erg.hoch);

    return erg;  // neu verrechnet
  }
}



// Hilfsklasse: speichert ein korresp. Paar (Zahl und Polynom), bzw. die Produkte dieser Paare, d.h. eine Liste von (a,b)
// Nach dem Verrechnungsschritt kann auch ein Produkt damit repräsentiert werden
class GNFS_Paar
{
  private GNFS_Teil_Ganzzahl erstens;  // a + b * m
  private GNFS_Teil_Ganzzahl zweitens;  // N(a + b * alpha)
  // die folgende Liste soll es später ermöglichen, das Produkt von (a + b * alpha) zusammenzustellen
  private List liste_a = null;
  private List liste_b = null;
  private final boolean DEBUG3 = false;  // Anzeige Polynomerstellung aus Z[alpha]

  // Konstruktor: setzt initial ein Basispaar (a + b * m) und (a + b * alpha)
  public GNFS_Paar (GNFS_Teil_Ganzzahl t1, GNFS_Teil_Ganzzahl t2, long a, long b)
  {
    erstens = t1; zweitens = t2;
    liste_a = new ArrayList(); liste_a.add(new Long(a));
    liste_b = new ArrayList(); liste_b.add(new Long(b));
  }

  // Kontruktor für Produkte, welche verrechnet wurden
  public GNFS_Paar (GNFS_Teil_Ganzzahl t1, GNFS_Teil_Ganzzahl t2, List a1, List b1, List a2, List b2)
  {
    erstens = t1; zweitens = t2;
    liste_a = new ArrayList(); liste_a.addAll(a1); liste_a.addAll(a2);
    liste_b = new ArrayList(); liste_b.addAll(b1); liste_b.addAll(b2);
    // Achtung: ist mit addAll auch gesagt, dass die Werte a und b korrespondierend eingefügt werden ?
  }


  // Liefert den höchsten belegten Index über der Faktorbasis zurück. Die von Term zweitens werden mit höherer Prior.
  // zurückgeliefert. Dies bedeutet, nach Aufruf muss festgestellt werden,
  // welcher Teil die höchste Variablenposition besitzt
  public int hoechstens()
  {
    for (int i = GNFS1.primBasis.length-1; i >= 0; --i)
    {
      // 1. Priorität zur Primzahl an i: die Norm
      if (zweitens.hoch.get(i)) return (2 * i + 1);  // aus Norm(a + b * alpha)
      if (erstens.hoch.get(i)) return (2 * i + 0);  // aus (a + b * m)
    }

    return -1;  // keine Exponenten gesetzt
  }


  // verrechnet die beiden Relationen paarweise
  public GNFS_Paar verrechnen(GNFS_Paar mit)
  {
    int w1 = this.hoechstens(); int w2 = mit.hoechstens();
    // kurze Prüfung
    if (w1 != w2)
    {
      // müssen beide die gleiche oberste Variable haben
      System.out.println ("Verrechnung Bew. " + w1 + " und " + w2 + ", obwohl nicht passend");
      System.exit(1);
    }
    if (w1 < 0 || w2 < 0) 
    {
      // wenn eine der beiden schon fertig ist, weil keine Variablen mehr enth.
      System.out.println("Verrechnung " + w1 + " und " + w2 + " nicht nötig");
      GNFS1.warten();
    }

    GNFS_Paar erg = new GNFS_Paar(this.erstens.verrechnen(mit.erstens), this.zweitens.verrechnen(mit.zweitens),
        this.liste_a, this.liste_b, mit.liste_a, mit.liste_b); 

    return erg;  // neu verrechnet
  }

  // macht Ausgabe
  public void ausgabe()
  {
    erstens.ausgabe(); 
    zweitens.ausgabe();
  }

  // liefert Wurzel von erstens, oder von zweitens bei algebr
  public long wurzel_liefern(boolean algebr)
  {
    return (algebr ? zweitens : erstens).wurz;
  }

  // liefert das aktuelle Polynom von zweitens. Rechnet modulo dem Hauptpolynom und n
  // Vorbedingung: Polynomgrad >= 2, wird bei Polynomerzeugung abgetestet
  public long [] polynom_liefern()
  {
    if ((liste_a.size() != liste_b.size()) || (liste_a.size() < 1)) { System.out.println ("Ungleiche Anzahl Elemente"); System.exit(1); }
    System.out.println ("Anzahl Paare zu verarbeiten: " + liste_a.size() );

    // Rechne modulo n, und modulo dem Vergleichspolynom, was gewählt wurde
    long n = GNFS1.holeN();
    long vergl [] = GNFS1.holePoly();
    int grad = GNFS1.holeGrad();
    if (vergl[grad] == 0) { System.out.println("Grad stimmt nicht"); System.exit(1); }

    long erg [] = new long[grad + 2];  // Grad 1, also 2 Koeff + 1 Reserve
    erg[0] = ((Long) liste_a.get(0)).longValue() % n; 
    erg[1] = ((Long) liste_b.get(0)).longValue() % n;  // a + b * alpha, d.h. es ergibt Z[alpha]
    int grad_aktuell = 1;
    if (DEBUG3) System.out.println ("Listen start a=" + erg[0] + " und b=" + erg[1] );
    for (int i = 1; (i < liste_a.size()) && (grad_aktuell >= 0); ++i)
    {
      // neues a
      long teila [] = erg.clone(); teila[grad_aktuell+1] = 0 /* zur Sicherh. */; 
      long aneu = ((Long) liste_a.get(i)).longValue();
      for (int j = grad_aktuell; j >= 0; --j) teila[j] = (erg[j] * aneu) % n;
      // neues b
      long teilb [] = erg.clone(); teilb[0] = 0;
      long bneu = ((Long) liste_b.get(i)).longValue();
      for (int j = grad_aktuell; j >= 0; --j) teilb[j+1] = (erg[j] * bneu) % n;  // +1, weil an b alpha hängt, damit verschoben

      // zusammenlegen, so dass es summe gibt, und noch modulo
      if (DEBUG3) System.out.println ("Neu a=" + aneu + " und b=" + bneu );
      ++grad_aktuell;
      for (int j = grad_aktuell; j >= 0; --j) erg[j] = (teila[j] + teilb[j]) % n;
      if (grad == grad_aktuell)  // wenn gl
      {
        //System.out.println ("Ergebnis:");
        //for (int k = grad; k >= 0; --k) System.out.print(erg[k] + " "); System.out.println("");
        //System.out.println ("Vergleich:");
        //for (int k = grad; k >= 0; --k) System.out.print(vergl[k] + " "); System.out.println("");

        // ein Vielfaches des festen Ausgangspolynoms (vergl) abziehen, modulo rechnen
        long inv = 1;
        try 
        {
          // wenn ggT(koeff, n) > 1, wirft modInverse eine Ausnahme
          inv = BigInteger.valueOf(vergl[grad]).modInverse(BigInteger.valueOf(n)).longValue();
        } 
        catch (Exception ex) 
        { 
          System.out.println ("ggT beim höchsten Koeff. " + BigInteger.valueOf(erg[grad]).gcd(BigInteger.valueOf(n)) ); System.exit(0);
        }
        long faktor = (erg[grad] * inv) % n;
        for (int j = grad; j >= 0; --j) erg[j] = ((erg[j] - faktor * vergl[j]) % n);
        if (erg[grad] != 0) { System.out.println ("modulo bei Polynom mit Fehler"); System.exit(1); }
        // Neuen aktuellen Grad ermitteln
        for (; grad_aktuell >= 0; --grad_aktuell) if (erg[grad_aktuell] != 0) break;
      } // Ende max. Grad erreicht
    }

    return erg;
  }

}