unit Polynome;

////////////////////////////////////////////////////////////////////////////////
//
//  Polynome
//  (c) 2007 Jakob Schöttl
//  Scho..ttl
//
//  TPolynom kapselt alle wichtigen Rechenoperationen
//  von mathematischen Polynomen; zusätzlich BruchPolynome.
//
////////////////////////////////////////////////////////////////////////////////

interface

uses Classes, SysUtils, Math, ExtCtrls, StdCtrls, Menus, Windows;

const
  maxGrad = 20;
  minGrad = 0;
  DefaultGrad = 3;

  //Error Messages (em)
  emInvalidPolynomGrad = 'Ungültiger PolynomGrad.';
  emInvalidGliederIndex = 'Ungültiger GliederIndex.';
  //function emParamException(...) //siehe weiter unten

  PolynomVariable = 'x';
  PolynomProduktOperator = '*';
  PolynomPotenzOperator = '^';
  PlusOperator = '+';
  MinusOperator = '-';

type

  TPolynomGlied = record
    Faktor: extended;
    Exponent: integer; //ReadOnly, sonst wäre es kein Polynom!
  end;

  TPolynomGlieder = array of TPolynomGlied; //index 0 bis FGliederCount -> rechts nach links

  TOutMathExprKind = (WithProduktOperator);
  TOutMathExprFlags = set of TOutMathExprKind;

  TInPolynomKind = (AutoSetGrad);
  TInPolynomFlags = set of TInPolynomKind;

  TNullstelle = record
    Wert: extended;
    Mehrfach: Boolean;
//    Vielfachheit: integer;
  end;
  TNullstellen = array of TNullstelle;

  TBruchPolynom = class;

  TPolynom = class(TObject)
  private
    FGrad: integer;
    FGlieder: TPolynomGlieder;
    FGliederCount: integer;
    FStoreOnChangingGrad: Boolean;
    FDefaultFaktor: extended;
    FInPolynomFlags: TInPolynomFlags;
    FOutMathExprFlags: TOutMathExprFlags;
  protected
    procedure SetGrad(Value: integer);
    function GetGlieder(index: integer): TPolynomGlied;
    procedure SetGlieder(index: integer; Value: TPolynomGlied);
    class procedure NsArrAufstSort(var n: TNullstellen); //Nullstellen-Array aufsteigend sortieren (Bubble-Sort)
  public
    constructor Create(aGrad: integer = DefaultGrad);
    destructor Destroy;
    function OutPolynom: string;
    function OutPolynomAsMathExpr: string;
    function OutPolynomForoooMath(With_newline: Boolean = True): string; //OpenOffice.org Math
    function OutPolynomForLATEX: string; //LATEX
    function InPolynom(PolynomStr: string; out Error: string): Boolean;
    function InPolynomAsMathExpr(PolynomStr: string; out Error: string): Boolean;

    property Glieder[index: integer]: TPolynomGlied read GetGlieder write SetGlieder;

    function ErgebnisFuer(x: extended): extended; 

    function ErsteAbleitung: TPolynom;
    function ZweiteAbleitung: TPolynom;
    function nteAbleitung(n: integer): TPolynom;
    function Stammfunktion: TPolynom;
    function Integral: TPolynom; overload;

    function Integral(UntereGrenze,ObereGrenze: extended): extended; overload;
    function FlaecheUnterKurve(UntereGrenze,ObereGrenze: extended): extended;

    function Nullstellen: TNullstellen;

//    {$DEFINE EXAKTELAENGE}
    {$UNDEF EXAKTELAENGE}
    {$IFDEF EXAKTELAENGE}
    function Laenge(UntereGrenze,ObereGrenze: extended): extended;     //nicht praktizierbar, weil ich die wurzel aus einem Polynom nicht ziehen kann
    {$ELSE}
    function Laenge(UntereGrenze,ObereGrenze: extended; AnzahlUnterteilungen: integer = 1000): extended;
    {$ENDIF}

    function Summe(Summand: TPolynom): TPolynom; overload;
    function Summe(Summand: extended): TPolynom; overload;
    function Differenz(Subtrahend: TPolynom): TPolynom; overload;
    function Differenz(Subtrahend: extended): TPolynom; overload;
    function Produkt(Faktor: TPolynom): TPolynom; overload;
    function Produkt(Faktor: extended): TPolynom; overload;
    function Quotient(Divisor: TPolynom; RestPolynom: TBruchPolynom): TPolynom; overload;
    function Quotient(Divisor: extended): TPolynom; overload;
    function Potenz(Exponent: integer): TPolynom; overload;
    function Negation: TPolynom; overload;


   { procedure Summe(Summand1,Summand2: TPolynom); overload;
    procedure Differenz(Minuend,Subtrahend: TPolynom); overload;
    procedure Produkt(Faktor1,Faktor2: TPolynom); overload;
    procedure Quotient(Dividend,Divisor: TPolynom); overload;
    procedure Potenz(Polynom: TPolynom; Exponent: integer); overload;   }

    function Equals(OtherPolynom: TPolynom): Boolean;
    function IsZero: Boolean;
    procedure Assign(Source: TPolynom);

    procedure FillFaktorenWith(f: extended);
    procedure Trim;

  published
    property Grad: integer read FGrad write SetGrad;
    property DefaultFaktor: extended read FDefaultFaktor write FDefaultFaktor;
    property StoreOnChangingGrad: Boolean read FStoreOnChangingGrad write FStoreOnChangingGrad;
    property InPolynomFlags: TInPolynomFlags read FInPolynomFlags write FInPolynomFlags;
    property OutMathExprFlags: TOutMathExprFlags read FOutMathExprFlags write FOutMathExprFlags;
  end;

  TBruchPolynom = class(TObject)
  private
    FZaehler: TPolynom;
    FNenner: TPolynom;
  public
    constructor Create;
    destructor Destroy; override;
    procedure Trim;

    function Nullstellen: TNullstellen;

    function Equals(OtherPolynom: TBruchPolynom): Boolean;
    function IsZero: Boolean;

    function ErsteAbleitung: TBruchPolynom;    //ungetestet!
//    function ZweiteAbleitung: TBruchPolynom;
//    function nteAbleitung(n: integer): TBruchPolynom;

    function OutPolynom: string;
    function OutPolynomAsMathExpr: string;
    function OutPolynomForoooMath(With_newline: Boolean = True): string; //OpenOffice.org Math
    function OutPolynomForLATEX: string; //LATEX
    property Zaehler: TPolynom read FZaehler;
    property Nenner: TPolynom read FNenner;
  end;


function PolynomGlied(Faktor: extended = 1): TPolynomGlied;
function ToPolynom(Faktoren: array of const): TPolynom;

implementation

function emParamException(ParamName,Meldung: string): string;
begin
  Result := 'Fehler in Parameter (' + ParamName + '): ' + Meldung;
end;







//Erstes Element in Liste wird zu Glied Nr. 1 bzw. Index 0 (Index somit gleich Exponent)!
function ToPolynom(Faktoren: array of const): TPolynom;
var i,l: integer;
begin
  Result := TPolynom.Create;
  l := Length(Faktoren) -1;
  Result.Grad := l;

  for i := 0 to l do
    With Faktoren[i] do
      case VType of
        vtExtended: Result.Glieder[i] := PolynomGlied(VExtended^);
        vtInteger: Result.Glieder[i] := PolynomGlied(VInteger);
        else begin
          FreeAndNil(Result);
          raise EMathError.Create(emParamException('Faktoren[' + IntToStr(i) + ']','Ungültiger ParameterType.'));
        end;
      end;
end;

function PolynomGlied(Faktor: extended = 1): TPolynomGlied;
begin
  Result.Exponent := -1;
  Result.Faktor := Faktor;
end;













function Seperate(Str: string; Seperator: char; out Part1,Part2: string): Boolean; overload;
var l,n: integer;
begin
  Result := False;
  l := Length(Str);
  for n := 1 to l do
    if Str[n] = Seperator then begin
      Part1 := Copy(Str,1,n -1);
      Part2 := Copy(Str,n +1,l - n);
      if (Part1 <> EmptyStr) and (Part2 <> EmptyStr) then Result := True;
      Break;
    end;
end;

type
  TStrParts = array of string;

function Seperate(Str: string; Seperator: char; out StrParts: TStrParts): Boolean; overload;
var n,l,i,j: integer;
    s: string;
begin
  Result := False;
  Str := Str + Seperator;
  l := Length(Str);     //Länge von Str
  j := 0;               //Soll-Länge von StrParts
  i := 0;               //Copy-Index
  for n := 1 to l do begin
    if Str[n] = Seperator then begin
      s := Copy(Str,i +1,n -i -1);
      i := n;
      if s <> EmptyStr then begin
        SetLength(StrParts,j +1);
        j := Length(StrParts);
        StrParts[j -1] := s;
      end;
    end;
  end;
  if j > 1 then Result := True;
end;

function Negativ(r: extended): Boolean; //True, wenn r < 0; False, wenn r >= 0;
begin
  if r < 0 then Result := True else Result := False;
end;

function _Round(const Wert: extended; const Stellen: integer): extended;
var
  Faktor,
  r: extended;
begin
  Faktor := IntPower(10, Stellen);

  if Wert > 0 then
    r := int(Wert * Faktor + 0.5)
  else
    r := int(Wert * Faktor - 0.5);

  r := r / Faktor;
  Result := r;
end;


















class procedure TPolynom.NsArrAufstSort(var n: TNullstellen);
var l,i,j: integer;
    r: extended;
begin
  //Nullstellen-Array aufsteigend sortieren (BubbleSort)
  l := Length(n);
  for i := 0 to l - 2 do
    for j := i +1 to l - 1 do
      if n[j].Wert < n[i].Wert then
      begin
        r := n[i].Wert;
        n[i].Wert := n[j].Wert;
        n[j].Wert := r;
      end;
end;

constructor TPolynom.Create(aGrad: integer = DefaultGrad);
begin
  FStoreOnChangingGrad := False;
  Grad := aGrad;
  FStoreOnChangingGrad := True;
  FDefaultFaktor := 1;
end;

destructor TPolynom.Destroy;
begin

  inherited;
end;


function TPolynom.OutPolynom: string;
var
  i: Integer;
  s: string;
begin
  Result := '';
  for i := 0 to FGliederCount - 1 do
    With FGlieder[i] do
    begin
      s := PolynomProduktOperator + PolynomVariable + PolynomPotenzOperator; //*x^
      s := FloatToStr(FGlieder[i].Faktor) + s;

      if (s[1] <> MinusOperator) and (i < Grad) then s := PlusOperator + s;
//      if Exponent < 0 then s := s + '(' + IntToStr(Exponent) + ')'           //wenn negativ einklammern
      { else }s := s + IntToStr(Exponent);
      Result := s + Result;
    end;
end;

function TPolynom.OutPolynomForLATEX: string; //LATEX
var n: integer;
begin
  Result := OutPolynomAsMathExpr;
  n := 0;
  While n < Length(Result) do begin
    if Result[n+1] = '*' then Delete(Result,n+1,1)
     else Inc(n);
  end;
end;

function TPolynom.InPolynom(PolynomStr: string; out Error: string): Boolean;
var sp: TStrParts;
    i,n,l: integer;
    p1,p2: string;
    pglieder: TPolynomGlieder;
    f: extended;
    e: integer;
    machen: Boolean;

begin
  Error := '';
  Result := True;
  SetLength(pglieder,FGliederCount);

  i := 2;
  While i < Length(PolynomStr) do
  begin
    if (PolynomStr[i] in [PlusOperator,MinusOperator]) then
    begin
      Insert(';',PolynomStr,i);
      Inc(i);
    end;
    Inc(i);
  end;

  Seperate(PolynomStr,';',sp); //egal ob true oder false, weil sonst fehler bei Grad=1

  machen := False;
  l := Length(sp);
  if l <> FGliederCount then
  begin
    if AutoSetGrad in FInPolynomFlags then
    begin
      Grad := l -1;
      machen := True;
    end else Error := 'Ungültige Anzahl der PolynomGlieder.';

  end else machen := True;



  if machen then
  begin
    n := FGliederCount;
    for i := 0 to FGliederCount - 1 do
    begin
      Dec(n);

      if not Seperate(StringReplace(sp[i],PolynomProduktOperator + PolynomVariable + PolynomPotenzOperator,';',[]),';',p1,p2) then
      begin
        Error := 'Fehler in PolynomGlied (Index ' + IntToStr(FGliederCount - i - 1) + ').';
        Break;
      end else
      begin
        if TryStrToFloat(p1,f) then pglieder[n] := PolynomGlied(f) else
        begin
          Error := 'Fehler in Faktor von PolynomGlied (Index ' + IntToStr(n) + ').';
          Break;
        end;
        if not (TryStrToInt(p2,e) and (FGlieder[n].Exponent = e)) then
        begin
          Error := 'Fehler in Exponent von PolynomGlied (Index ' + IntToStr(n) + ').';
          Break;
        end;
      end;
    end;

  end;

  if Error <> '' then Result := False;

  if Result then
    for i := 0 to FGliederCount - 1 do
    begin
      Glieder[i] := PolynomGlied(pglieder[i].Faktor);
    end;

end;

function TPolynom.InPolynomAsMathExpr(PolynomStr: string; out Error: string): Boolean;
var i,l,n: integer;
begin
  l := Length(PolynomStr);
  i := 1;
  While i <= l do
  begin
    case PolynomStr[i] of
      ' ':
        begin
          Delete(PolynomStr,i,1);
          Dec(i);
          Dec(l);
        end;
      PolynomVariable:
        if (i <= 1) or (PolynomStr[i -1] <> PolynomProduktOperator) then
        begin
          Insert(PolynomProduktOperator,PolynomStr,i);
          Dec(i);
          Inc(l);
        end else
        if (i < l) and (PolynomStr[i +1] <> '^') then
        begin

          Insert('^1',PolynomStr,i+1);
          Inc(l,2);
        end;
      PolynomProduktOperator:
        if (i <= 1) or not (PolynomStr[i -1] in ['1','2','3','4','5','6','7','8','9','0',',','E','e']) then
        begin
          Insert('1',PolynomStr,i);
          Inc(l);
        end;
    end;
    Inc(i);
  end;

  l := Length(PolynomStr);
  for i := l downto 1 do
    if PolynomStr[i] = '^' then Break;
  if not TryStrToInt(Copy(PolynomStr,i+1,l),n) then
    PolynomStr := PolynomStr + '*x^0';


  Result := InPolynom(PolynomStr,Error);
end;

procedure TPolynom.SetGrad(Value: integer);
var
  i: Integer;
  pgr: TPolynomGlieder;
begin
  if Value <> FGrad then
    if (Value <= maxGrad) and (Value >= minGrad) then
    begin

      if FStoreOnChangingGrad then
        if Value > FGrad then
        begin
          SetLength(pgr,FGrad +1);
          for i := 0 to FGrad do
            pgr[i] := Glieder[i];
        end
         else
         begin
           SetLength(pgr,Value +1);
           for i := 0 to Value do
             pgr[i] := Glieder[i];
         end;


      FGrad := Value;
      FGliederCount := FGrad + 1;
      SetLength(FGlieder,FGliederCount);

      for i := 0 to FGliederCount - 1 do
        With FGlieder[i] do
        begin
          Faktor := FDefaultFaktor;
          Exponent := i;
        end;


      if FStoreOnChangingGrad then
        for i := Low(pgr) to High(pgr) do
          FGlieder[i] := pgr[i];


    end else raise EMathError.Create(emInvalidPolynomGrad);
end;

function TPolynom.GetGlieder(index: integer): TPolynomGlied;
begin
  if (index < FGliederCount) and (index >= 0) then
  begin
    Result := FGlieder[index];
  end else raise EAccessViolation.Create(emInvalidGliederIndex);
end;

procedure TPolynom.SetGlieder(index: integer; Value: TPolynomGlied);
begin
  if (index < FGliederCount) and (index >= 0) then
  begin
    FGlieder[index].Faktor := Value.Faktor;
//    FGlieder[index].Exponent := Value.Exponent;
  end else raise EAccessViolation.Create(emInvalidGliederIndex);
end;

function TPolynom.ErgebnisFuer(x: extended): extended;
var i: integer;
begin
  Result := 0;
  for i := 0 to FGliederCount - 1 do
    With Glieder[i] do
      Result := Result + (Faktor * Power(x,Exponent));

end;

function TPolynom.Integral(UntereGrenze,ObereGrenze: extended): extended;
{$DEFINE RICHTIGESINTEGRAL}    //nicht-"richtiges integral" wäre die wirkliche fläche zwischen dem Graphen einer Funktion und der x-Achse
{$IFDEF RICHTIGESINTEGRAL}
begin
  With Integral do
  try
    Result := ErgebnisFuer(ObereGrenze) - ErgebnisFuer(UntereGrenze);
  finally
    Free;
  end;
end;
{$ELSE}    //folgendes wäre nur die fläche zwischen x-Achse und Kurve!
var p: TPolynom;
    n: TNullstellen;
    a: array of extended;
    i,j,l: integer;
    r: extended;
    merker: extended;
begin

  if UntereGrenze > ObereGrenze then
  begin
    //Obere- und UntereGrenze vertauschen (für korrekte Berechnung)
    merker := UntereGrenze;
    UntereGrenze := ObereGrenze;
    ObereGrenze := merker;

  end else
  if UntereGrenze = ObereGrenze then  //Sind Grenzwerte gleich? -> Fläche = 0
  begin
    Result := 0;
    Exit;
  end;


  p := Integral;

  //Die Fläche unter dem Graph ist die Summe aller Beträge der Flächen, angefangen bei UnterGrenze über alle dazwischenliegenden Nullstellen (aufsteigend sortiert) hin bis zu ObereGrenze.
  //-------------------------------------------------------------------------------------------------------------------------------------------------------------------

  //Nullstellen berechnen
  n := Nullstellen;

  NsArrAufstSort(n);

  //array a für Berechnungen:

  //Untere Grenze zum arr hinzufügen
  SetLength(a,1);
  a[0] := UntereGrenze;

  //dazwischenliegende Nullstellen zwischen Unterer und Oberer Grenze einfügen
  for i := 0 to l - 1 do //l noch unverändert von oben (länge des Nullstellen-Arrays)
    if (n[i].Wert > UntereGrenze) and (n[i].Wert < ObereGrenze) then
    begin
      l := Length(a);
      SetLength(a,l + 1);
      a[l] := n[i].Wert;
    end else Break;

  //Obere Grenze zum arr hinzufügen
  l := Length(a);
  SetLength(a,l + 1);
  a[l] := ObereGrenze;

  //Summe der Beträge der Flächen berechnen
  Result := 0;
  for i := 0 to l - 1 do //l = länge von a - 1; nach algo a.count -2, also hier: ... -1
    Result := Result + abs(p.ErgebnisFuer(a[i+1]) - p.ErgebnisFuer(a[i]));

//  Result := p.ErgebnisFuer(ObereGrenze) - p.ErgebnisFuer(UntereGrenze);
  p.Free;

end;
{$ENDIF}



function TPolynom.FlaecheUnterKurve(UntereGrenze,ObereGrenze: extended): extended;
var p: TPolynom;
    n: TNullstellen;
    a: array of extended;
    i,j,l: integer;
    r: extended;
    merker: extended;

begin

  if UntereGrenze > ObereGrenze then
  begin
    //Obere- und UntereGrenze vertauschen (für korrekte Berechnung)
    merker := UntereGrenze;
    UntereGrenze := ObereGrenze;
    ObereGrenze := merker;

  end else
  if UntereGrenze = ObereGrenze then  //Sind Grenzwerte gleich? -> Fläche = 0
  begin
    Result := 0;
    Exit;
  end;


  p := Integral;

  //Die Fläche unter dem Graph ist die Summe aller Beträge der Flächen, angefangen bei UnterGrenze über alle dazwischenliegenden Nullstellen (aufsteigend sortiert) hin bis zu ObereGrenze.
  //---------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------

  //Nullstellen berechnen
  n := Nullstellen;

  //array a für Berechnungen:

  //Untere Grenze zum arr hinzufügen
  SetLength(a,1);
   a[0] := UntereGrenze;
   l := Length(n);

   //dazwischenliegende Nullstellen zwischen Unterer und Oberer Grenze einfügen
   for i := 0 to l - 1 do
     if (n[i].Wert > UntereGrenze) and (n[i].Wert < ObereGrenze) then
     begin
       l := Length(a);
       SetLength(a,l + 1);
       a[l] := n[i].Wert;
     end;

   //Obere Grenze zum arr hinzufügen
   l := Length(a);
   SetLength(a,l + 1);
   a[l] := ObereGrenze;

   //Summe der Beträge der Flächen berechnen
   Result := 0;
   for i := 0 to l - 1 do //l = länge von a - 1; nach algo a.count -2, also hier: ... -1
     Result := Result + abs(p.ErgebnisFuer(a[i+1]) - p.ErgebnisFuer(a[i]));

   //  Result := p.ErgebnisFuer(ObereGrenze) - p.ErgebnisFuer(UntereGrenze);
   p.Free;
   
end;
















function TPolynom.Nullstellen: TNullstellen;

  function NullstellenFinden(Polynom: TPolynom): TNullstellen;

    procedure AddNullstelle(var Nullstellen: TNullstellen; x: extended);
    var l: integer;
    begin
      l := Length(Nullstellen);
      SetLength(Nullstellen,l + 1);
      Nullstellen[l].Wert := x;
    end;

  const
    Naeherung_msTimeOut = 2000;
    Naeherung_StellenNachKomma = 10;
    Naeherung_epsilon = 1E-10;
    //E-10 entspricht 10 Stellen nach Komma, also Runden auf 10 Stellen nach Komma

    function Naeherungswert_TangentenVerfahren(x: extended): extended;
    var
      a: integer;
      x_neu: extended;
      Ableitung: TPolynom;

    begin
      a := GetTickCount;

      x_neu := x;

      Ableitung := Polynom.ErsteAbleitung;

      repeat

        x := x_neu;
        x_neu := x - Polynom.ErgebnisFuer(x) / Ableitung.ErgebnisFuer(x);

      until (x_neu = x) or (abs(x_neu - x) < Naeherung_epsilon) {or ((GetTickCount - a) > Naeherung_msTimeOut)};

      Ableitung.Free;

      Result := x_neu;
    end;

    function Naeherungswert_SehnenVerfahren(x_1,x_2: extended): extended;
    var
      a: integer;
      x,x_neu: extended;
      fvonx_2,fvonx_neu: extended;

    begin
      a := GetTickCount;

      repeat

        x_neu := x_1 - Polynom.ErgebnisFuer(x_1) *( (x_2 - x_1) /
                                                    (Polynom.ErgebnisFuer(x_2) - Polynom.ErgebnisFuer(x_1)) );


        fvonx_neu := Polynom.ErgebnisFuer(x_neu);

        if fvonx_neu = 0 then Break else

        //Haben Funktionswerte gleiches Vorzeichen?
        if not (Negativ(Polynom.ErgebnisFuer(x_1)) xor Negativ(fvonx_neu)) then
        begin
          x := x_1;
          x_1 := x_neu;
        end else
//        if not (Negativ(Polynom.ErgebnisFuer(x_2)) xor Negativ(fvonx_neu)) then //if-Anweisung unnötig, weil's gar nicht anders sein kann.
        begin
          x := x_2;
          x_2 := x_neu;
        end;


      until (x_neu = x) or (abs(x_neu - x) < Naeherung_epsilon) {or ((GetTickCount - a) > Naeherung_msTimeOut)};

      Result := x_neu;
    end;



  var
    x,Radiant: extended; //nur für mitternachtsformel
    y_1,y_2: extended;
    x_1,x_2: extended;
    Ableitung: TPolynom;
    N: TNullstellen;
    sIndex,eIndex: integer; //start-/end-Index für HP-/TP-Array
    l,i: integer;
  begin

    SetLength(Result,0);

    case Polynom.Grad of
      0: {nichts} ;
      1: //Lineare Polynome: y=mx+t <=> x_0=-t/m
         AddNullstelle(Result,- Polynom.Glieder[0].Faktor / Polynom.Glieder[1].Faktor);
      2: //Quadratische Polynome: Mitternachtsformel
         begin

           Radiant := Power(Polynom.Glieder[1].Faktor,2) - 4*Polynom.Glieder[2].Faktor*Polynom.Glieder[0].Faktor;
           if Radiant >= 0 then
           begin
             AddNullstelle(Result,(- Polynom.Glieder[1].Faktor + sqrt(Radiant)) / (2*Polynom.Glieder[2].Faktor));

             //Radiant hat sich nicht geändert
             x := (- Polynom.Glieder[1].Faktor - sqrt(Radiant)) / (2*Polynom.Glieder[2].Faktor);
             if x <> Result[0].Wert then
               AddNullstelle(Result,x);
           end;
         end;

      else //Polynome mit Grad > 2
      begin


        //Durch rekursives Nullstellen-Finden der Ableitung (...) Hoch- und Tiefpunkte finden
        Ableitung := Polynom.ErsteAbleitung;
        N := NullstellenFinden(Ableitung);
        Ableitung.Free;
        NsArrAufstSort(N);


        sIndex := 0;
        eIndex := Length(N) -1;

        //links: Näherung mit TangentenVerfahren
        y_1 := Polynom.ErgebnisFuer(N[0].Wert);
        y_2 := Polynom.ErgebnisFuer(N[0].Wert - 1);
        if _Round(y_1,Naeherung_StellenNachKomma) = 0 then
        begin
          AddNullstelle(Result,N[0].Wert);
          Inc(sIndex);
        end else
        if ((y_2 < y_1) and (y_1 > 0)) or
           ((y_2 > y_1) and (y_1 < 0)) then
             AddNullstelle(Result,_Round(Naeherungswert_TangentenVerfahren(N[0].Wert - 1),Naeherung_StellenNachKomma));


        //rechts: Näherung mit TangentenVerfahren
        y_1 := Polynom.ErgebnisFuer(N[eIndex].Wert);
        y_2 := Polynom.ErgebnisFuer(N[eIndex].Wert + 1);
        if _Round(y_1,Naeherung_StellenNachKomma) = 0 then
        begin
          AddNullstelle(Result,N[eIndex].Wert);
          Dec(eIndex);
        end else
        if ((y_2 < y_1) and (y_1 > 0)) or
           ((y_2 > y_1) and (y_1 < 0)) then
             AddNullstelle(Result,_Round(Naeherungswert_TangentenVerfahren(N[eIndex].Wert + 1),Naeherung_StellenNachKomma));



        //mitte: Näherung mit SehnenVerfahren
        i := sIndex;
        While i < eIndex do
        begin
          x_1 := N[i].Wert;
          x_2 := N[i + 1].Wert;

          y_1 := Polynom.ErgebnisFuer(x_1);
          y_2 := Polynom.ErgebnisFuer(x_2);

  //        if y_1 = 0 then AddNullstelle(Result,x_1) else //Dieser Fall sollte nicht eintreten können! Weil x_1 schon bei der näherung links ausgeschlossen wurde. Und zwei HP/TP nebeneinander, für die beide gilt y=0 gibt es nicht.
          if _Round(y_2,Naeherung_StellenNachKomma) = 0 then
          begin
            AddNullstelle(Result,x_2);
            Inc(i);
          end else

              //eins negativ, das andere positiv
          if Negativ(y_1) = not Negativ(y_2) then //unterschiedliche Vorzeichen (0 kann theoretisch keine der beiden Zahlen sein)
            AddNullstelle(Result,_Round(Naeherungswert_SehnenVerfahren(N[i].Wert,N[i + 1].Wert),Naeherung_StellenNachKomma));


          







          Inc(i);  //inc schleifenvar
        end;

      end; //Ende von Grad > 2
    end; //Ende von Case
  end;












var
  Ableitung: TPolynom;
  n: integer;
  
begin
  Trim;
  Result := NullstellenFinden(Self);

  //Nullstellen sortieren (NsArrAufstSort)
  NsArrAufstSort(Result);

  //Nullstellen (Result) auf Mehrfachheit prüfen: f'(x_null) = 0
  Ableitung := ErsteAbleitung;
  for n := 0 to Length(Result) - 1 do
    With Result[n] do
      if Ableitung.ErgebnisFuer(Wert) = 0 then
        Mehrfach := True else
        Mehrfach := False;
  Ableitung.Free;
    
    
end;























{$IFDEF EXAKTELAENGE}
function TPolynom.Laenge(UntereGrenze,ObereGrenze: extended): extended;
begin
  //Länge einer Funktion f(x) im Intervall [a;b]:
  //Integral a bis b von (der Wurzel aus (der Summe 1 und (der abgeleitete Funktion zum Quadrat)))
  // oder
  //Integral a bis b von sqrt(1 + (f'(x)^2))
{  With Ableitung do
  try
    With Potenz(2) do
    try
      With Summe(1) do
      try
        With Wurzel(2) do  //hier ist das problem
        try
          Result := Integral(UntereGrenze,ObereGrenze);
        finally
          Free;
        end;
      finally
        Free;
      end;
    finally
      Free;
    end;
  finally
    Free;
  end; }
end;

{$ELSE}
function TPolynom.Laenge(UntereGrenze,ObereGrenze: extended; AnzahlUnterteilungen: integer = 1000): extended;
var diff,step,merker: extended;
    n: integer;
    lnegativ: Boolean;
begin

  if AnzahlUnterteilungen <= 0 then raise EMathError.Create(emParamException('AnzahlUnterteilungen','Anzahl der Unterteilungen muss größer sein als 0.'));


  //Ist UntereGrenze > ObereGrenze? -> Laenge ist negativ
  lnegativ := False;
  if UntereGrenze > ObereGrenze then
  begin
    lnegativ := True;  //wird am Ende der Funktion auf Laenge (Result) übertragen

    //Obere- und UntereGrenze vertauschen (für korrekte Berechnung)
    merker := UntereGrenze;
    UntereGrenze := ObereGrenze;
    ObereGrenze := merker;

  end else
  if UntereGrenze = ObereGrenze then  //Sind Grenzwerte gleich? -> Fläche = 0
  begin
    Result := 0;
    Exit;
  end;

  //Länge des Graphen einer Funktion in einem Intervall := Summe aller euklidischen Längen zwischen den dazwischenliegenden Punkten


  diff := ObereGrenze - UntereGrenze;
  step := diff / AnzahlUnterteilungen;
  Result := 0;

  for n := 0 to AnzahlUnterteilungen -1 do
    Result := Result +
      sqrt(Power(step,2) +
           Power(ErgebnisFuer(step * (n+1) + UntereGrenze) - ErgebnisFuer(step * n + UntereGrenze) ,2));


  if lnegativ then Result := -Result;

end;
{$ENDIF}

function TPolynom.ErsteAbleitung: TPolynom;
var i: integer;
    pg: TPolynomGlied;
begin
  Result := TPolynom.Create;

  if Grad > 0 then
  begin
    Result.Grad := Grad -1;

    for i := 1 to FGliederCount - 1 do
      With Glieder[i] do
      begin
        pg.Faktor := Exponent * Faktor;
        pg.Exponent := Exponent - 1;
        Result.Glieder[i -1] := pg;
      end;

  end else                   //Konstanten ableiten -> 0
   begin
     Result.Grad := 0;
     pg := Result.Glieder[0];
     pg.Faktor := 0;
     Result.Glieder[0] := pg;
   end;
end;

function TPolynom.ZweiteAbleitung: TPolynom;
begin
  Result := nteAbleitung(2);
end;

function TPolynom.Stammfunktion: TPolynom;
var i: integer;
    pg: TPolynomGlied;
begin
  Result := TPolynom.Create;

  Result.Grad := Grad +1;

  for i := 0 to FGliederCount - 1 do
    With Glieder[i] do
    begin
      pg.Exponent := Exponent + 1;
      pg.Faktor := Faktor / pg.Exponent;
      Result.Glieder[i +1] := pg;
    end;

  pg.Exponent := 0;
  pg.Faktor := 0;
  Result.Glieder[0] := pg;
end;

function TPolynom.Integral: TPolynom;
begin
  Result := Stammfunktion;
end;

function TPolynom.nteAbleitung(n: integer): TPolynom;    //ohne instanzen-array wäre viel schöner!
{$DEFINE OPTIMIERTER_ALGOLRITHMUS}
{$IFNDEF OPTIMIERTER_ALGOLRITHMUS} //If Not Defined
var
  i: Integer;
  parr: array of TPolynom;
begin
  Result := nil;
  if (n < 1) then raise EAccessViolation.Create(emAbleitungsFehler)
   else
   begin
     SetLength(parr,n + 1);
     parr[0] := Self;
     for i := 1 to n do
       parr[i] := parr[i - 1].ErsteAbleitung;

     for i := 1 to n - 1 do
       parr[i].Free;
     Result := parr[n];
     parr[n] := nil;
     parr[0] := nil;
   end;
end;
{$ELSE}
var
  i,j,x,SExp: integer;
  f: extended;
begin
  Result := TPolynom.Create;

  if n > Self.Grad then
    With Result do
    begin
      Grad := 0;
      Glieder[0] := Polynomglied(0);
    end else


    begin
      Result.Grad := Self.Grad - n;
      j := Self.Grad - Result.Grad;
      for i := 0 to Result.Grad do
      begin
        f := Self.Glieder[j].Faktor;
        SExp := Self.Glieder[j].Exponent;
        for x := SExp downto (SExp - n + 1) do
          f := f * x;

        Result.Glieder[i] := PolynomGlied(f);

        Inc(j);
      end;
    end;

end;
{$ENDIF}

function TPolynom.OutPolynomAsMathExpr: string;
var
  i: Integer;
  s,sFakt: string;
  l: integer;
begin
  Result := '';
  for i := 0 to FGliederCount - 1 do
    With FGlieder[i] do
      if Faktor <> 0 then
      begin


        if Exponent = 0 then s := '' else
        if Exponent = 1 then s := PolynomVariable { x }
        else
        begin
          s := PolynomVariable + PolynomPotenzOperator; // x^
//        if Exponent < 0 then s := s + '(' + IntToStr(Exponent) + ')'           //wenn negativ einklammern
          { else }s := s + IntToStr(Exponent);
        end;

        sFakt := FloatToStr(Faktor);
        if (abs(Faktor) <> 1) then
        begin
          if s <> '' then
          begin
            if WithProduktOperator in FOutMathExprFlags then s := sFakt + PolynomProduktOperator + s
             else s := sFakt + s;
          end
           else s := sFakt;
        end else
          if i = 0 then s := sFakt else
            if Faktor < 0 then s := MinusOperator + s;;

        if (s[1] <> MinusOperator) and (i < Grad) then s := PlusOperator + s;

        Result := s + Result;
      end;
  if Result = '' then Result := '0' else
  if Result[1] = PlusOperator then Delete(Result,1,1);
end;

function TPolynom.OutPolynomForoooMath(With_newline: Boolean = True): string;
var n: integer;
begin
  Result := OutPolynomAsMathExpr;
  n := 0;
  While n < Length(Result) do begin
    if Result[n+1] = '*' then Delete(Result,n+1,1)
     else Inc(n);
  end;
  if With_newline then Result := Result + ' newline';
end;

function TPolynom.Equals(OtherPolynom: TPolynom): Boolean;
var n: integer;
begin
  Result := True;
  if OtherPolynom.Grad = Grad then
  begin
    for n := 0 to FGliederCount -1 do
    begin
      if (OtherPolynom.Glieder[n].Faktor <> Glieder[n].Faktor)
//        or (OtherPolynom.Glieder[n].Exponent <> Glieder[n].Exponent)
         then
           Result := False;
    end;
  end else Result := False;
end;

function TPolynom.Summe(Summand: TPolynom): TPolynom;
var n: integer;
    pg: TPolynomGlied; //Temporärer Speicher
//    Summand1: TPolynom; //Summand; der Summand mit dem höheren Grad, zu dem Summand2 addiert wird
    Summand2: TPolynom; //Summand, der zum Ergebnis addiert wird

begin
  //Result := Result + Summand2
  Result := TPolynom.Create;
  if Self.Grad > Summand.Grad then begin
    Result.Assign(Self);
    Summand2 := Summand;
  end
   else
   begin
     Result.Assign(Summand);
     Summand2 := Self;
   end;

  for n := 0 to Summand2.Grad do
  begin
    pg := Result.Glieder[n];

    pg.Faktor := pg.Faktor + Summand2.Glieder[n].Faktor;
    Result.Glieder[n] := pg;
  end;
  Result.Trim;

end;

function TPolynom.Summe(Summand: extended): TPolynom;
begin
  Result := TPolynom.Create;
  Result.Assign(Self);
  Result.FGlieder[0].Faktor := Result.FGlieder[0].Faktor + Summand;  
end;

function TPolynom.Differenz(Subtrahend: TPolynom): TPolynom;
var n: integer;
    pg: TPolynomGlied; //Temporärer Speicher
//    Minuend: TPolynom; //Minuend mit dem höheren Grad, von dem Subtrahend2 abgezogen wird

begin
  //Result := Self;
  //Result := Result - Subtrahend;


  Result := TPolynom.Create;
  if Self.Grad >= Subtrahend.Grad then
    Result.Grad := Self.Grad
   else
     Result.Grad := Subtrahend.Grad;           //Grad von Result auf höheren Grad setzen

  Result.FillFaktorenWith(0);
  for n := 0 to Self.Grad do Result.Glieder[n] := PolynomGlied(Self.Glieder[n].Faktor); //Result := Self;

  for n := 0 to Subtrahend.Grad do             //Result := Result - Subtrahend;
    Result.Glieder[n] := PolynomGlied(Result.Glieder[n].Faktor - Subtrahend.Glieder[n].Faktor);

  Result.Trim;
end;

function TPolynom.Differenz(Subtrahend: extended): TPolynom;
begin
  Result := TPolynom.Create;
  Result.Assign(Self);
  Result.FGlieder[0].Faktor := Result.FGlieder[0].Faktor - Subtrahend;
end;

procedure TPolynom.Assign(Source: TPolynom);
var n: integer;
begin
  Grad := Source.Grad;
  for n := 0 to Grad do
    With Self.FGlieder[n] do
    begin
      Faktor := Source.Glieder[n].Faktor;
//      Exponent := Source.Glieder[n].Exponent;
    end;
end;

function TPolynom.Produkt(Faktor: TPolynom): TPolynom;
var h,i,j: integer;
begin
  Trim;
  Faktor.Trim;

  Result := TPolynom.Create;
  Result.Grad := Self.Grad + Faktor.Grad;

  for i := 0 to Result.Grad do
    Result.Glieder[i] := PolynomGlied(0);

  for i := 0 to Self.Grad do
    for j := 0 to Faktor.Grad do
    begin
      h := Self.Glieder[i].Exponent + Faktor.Glieder[j].Exponent;
      Result.Glieder[h] := PolynomGlied(Result.Glieder[h].Faktor + (Self.Glieder[i].Faktor) * Faktor.Glieder[j].Faktor);
    end;

  Result.Trim;
end;

function TPolynom.Produkt(Faktor: extended): TPolynom;
var n: integer;
begin
  Result := TPolynom.Create;
  Result.Assign(Self);
  for n := 0 to Result.Grad do
    Result.FGlieder[n].Faktor := Result.FGlieder[n].Faktor * Faktor;
end;

function TPolynom.Quotient(Divisor: TPolynom; RestPolynom: TBruchPolynom): TPolynom;
var P,            //(Minuend-)Polynom
    S,            //Subtrahend-Polynom
    F,            //Faktor-Polynom
    M: TPolynom;  //Merker (zum nachträglichen Freigeben)

    GDiff,             //Grad-Differez
    n,                 //Grad von Result (=GDiff)
    PDGrad,            //Divisor.Grad
    i: integer;        //Schleifenvariable

    PDFaktor: extended;    //Divisor.Glieder[Grad].Faktor

begin

  Divisor.Trim;
  if not Divisor.IsZero then //es wird nicht durch 0 geteilt
  begin

    Self.Trim;

    GDiff := Self.Grad - Divisor.Grad;
    if GDiff >= 0 then
    try
      Result := TPolynom.Create;
      Result.Grad := GDiff;
      Result.FillFaktorenWith(0);

      PDGrad := Divisor.Grad;
      PDFaktor := Divisor.Glieder[PDGrad].Faktor; //Faktor durch den immer geteilt wird (der mit dem höchsten Exponenten)

      P := TPolynom.Create;
      P.Assign(Self);

      F := TPolynom.Create;

//      n := GDiff;
      While P.Grad >= Divisor.Grad do
      begin

        n := P.Grad - PDGrad;

        Result.Glieder[n] := PolynomGlied(P.Glieder[P.Grad].Faktor / PDFaktor);

        F.Grad := n;
        for i := 0 to n - 1 do
          F.Glieder[i] := PolynomGlied(0);
        F.Glieder[n] := PolynomGlied(Result.Glieder[n].Faktor);

        S := F.Produkt(Divisor);

        M := P;
        P := P.Differenz(S);
        M.Free;
        S.Free;

//        Dec(n);
      end;

    finally

      RestPolynom.Zaehler.Assign(P);
      RestPolynom.Nenner.Assign(Divisor);

      P.Free;
      F.Free;

    end else
    begin
      Result := ToPolynom([0]);
      RestPolynom.Zaehler.Assign(Self);
      RestPolynom.Nenner.Assign(Divisor);
    end;


  end else
  begin
    raise EZeroDivide.Create('Divisor ist 0.');
    FreeAndNil(Result);
  end;
end;

function TPolynom.Quotient(Divisor: extended): TPolynom;
var n: integer;
begin
  if Divisor = 0 then
  begin
    raise EZeroDivide.Create('Divisor ist 0.');
    FreeAndNil(Result);
    Exit;
  end;

  Result := TPolynom.Create;
  Result.Assign(Self);
  for n := 0 to Result.Grad do
    Result.FGlieder[n].Faktor := Result.FGlieder[n].Faktor / Divisor;
end;

function TPolynom.Potenz(Exponent: integer): TPolynom;
var n: integer;
    m: TPolynom;
begin
  if Exponent < 0 then raise EMathError(emParamException('Exponent','Der Exponent darf nicht kleiner als 0 sein.')) else
  if Exponent = 0 then Result := ToPolynom([1]) else
  begin
//    if (Exponent * Grad) > maxGrad then raise EMathError.Create('Polynom.Grad-Overflow: Der Grad des Ergebnis-Polynom ist zu groß.');

    Result := ToPolynom([1]);
    for n := 0 to Exponent - 1 do
    begin
      m := Result;
      Result := Result.Produkt(Self);
      FreeAndNil(m);
    end;
  end;
end;

function TPolynom.Negation: TPolynom;
begin
  Result := Produkt(-1);
end;

function TPolynom.IsZero: Boolean;
var i: integer;
begin
  Result := True;
  for i := 0 to Grad do
    if Glieder[i].Faktor <> 0 then
    begin
      Result := False;
      Break;
    end;
end;

procedure TPolynom.Trim;      //schneidet die linken summanden (glieder) ab, wenn ihr Faktor 0 ist
var i,n: integer;
    tp: TPolynom; //Temporäres Polynom
begin
  if Grad > 0 then
  begin
    for i := Grad downto 1 do
      if Glieder[i].Faktor <> 0 then Break;

    if i < Grad then
    try
      tp := TPolynom.Create;
      tp.Grad := i;
      for n := 0 to tp.Grad do
        tp.Glieder[n] := PolynomGlied(Glieder[n].Faktor);
      Self.Assign(tp);
    finally
      tp.Free;
    end;

  end;
end;

procedure TPolynom.FillFaktorenWith(f: extended);
var i: integer;
begin
  for i := 0 to Grad do
    Glieder[i] := PolynomGlied(f);

end;


















// TBruchPolynom ///////////////////////////////////////////////////////

constructor TBruchPolynom.Create;
begin
  FZaehler := TPolynom.Create;
  With FZaehler do
  begin
    Grad := 0;
  end;
  FNenner := TPolynom.Create;
  With FNenner do
  begin
    Grad := 0;
  end;
end;

destructor TBruchPolynom.Destroy;
begin
  FZaehler.Free;
  FNenner.Free;
  inherited Destroy;
end;


function TBruchPolynom.ErsteAbleitung: TBruchPolynom;
var
  v,u,             //v und u (Zähler und Nenner)
  v_,u_: TPolynom; //v' und u'
  p1,p2: TPolynom;
begin
//  f(x)=u(x)/v(x) => f'(x)=(u'(x)*v(x)-u(x)*v'(x))/(v(x))^2
  try
    v := Zaehler;
    u := Nenner;
    v_ := Zaehler.ErsteAbleitung;
    u_ := Nenner.ErsteAbleitung;

    Result := TBruchPolynom.Create;

    //Zähler
    p1 := u_.Produkt(v);
    p2 := u.Produkt(v_);

    Result.FZaehler.Free;
    Result.FZaehler := p1.Differenz(p2);

    //Nenner
    Result.FNenner.Free;
    Result.FNenner := v.Potenz(2);

  finally
    v.Free;
    u.Free;
    v_.Free;
    u_.Free;
    p1.Free;
    p2.Free;
  end;
end;


procedure TBruchPolynom.Trim;
begin
  FZaehler.Trim;
  FNenner.Trim;
end;

function TBruchPolynom.OutPolynom: string;
begin
  Result := '(' + FZaehler.OutPolynom + ')/(' + FNenner.OutPolynom + ')';
end;

function TBruchPolynom.OutPolynomAsMathExpr: string;
begin
  Result := '(' + FZaehler.OutPolynomAsMathExpr + ')/(' + FNenner.OutPolynomAsMathExpr + ')';
end;

function TBruchPolynom.OutPolynomForoooMath(With_newline: Boolean = True): string; //OpenOffice.org Math
begin
  Result := '{' + FZaehler.OutPolynomForoooMath(False) + '} over {' + FNenner.OutPolynomForoooMath(False) + '}';
  if With_newline then Result := Result + ' newline';

end;

function TBruchPolynom.OutPolynomForLATEX: string; //LATEX
var n: integer;
begin
  Result := '\frac{' + Zaehler.OutPolynomAsMathExpr + '}{' + Nenner.OutPolynomAsMathExpr + '}';
  n := 0;
  While n < Length(Result) do begin
    if Result[n+1] = '*' then Delete(Result,n+1,1)
     else Inc(n);
  end;
end;

function TBruchPolynom.Nullstellen: TNullstellen;
begin
  Result := FZaehler.Nullstellen;
end;

function TBruchPolynom.Equals(OtherPolynom: TBruchPolynom): Boolean;
begin
  if FZaehler.Equals(OtherPolynom.FZaehler) and FNenner.Equals(OtherPolynom.FNenner) then Result := True
   else Result := False;
end;

function TBruchPolynom.IsZero: Boolean;
begin
  if FZaehler.IsZero {and not FNenner.IsZero} then Result := True
   else Result := False;  
end;

end.

