Einzelnen Beitrag anzeigen

Benutzerbild von FAlter
FAlter

Registriert seit: 21. Jul 2004
Ort: Ostfildern
1.096 Beiträge
 
FreePascal / Lazarus
 
#9

Re: Guter Algorithmus zum multiplizieren großer Zahlen gesuc

  Alt 15. Nov 2008, 17:38
Hi,

hier mein aktueller Stand, ich habe in Assembler eine Proc zum unsigned multiplizieren geschrieben, die sogar ein 256 bittiges Ergebnis zurückliefert. Und dann In Delphi die Behandlung des Vorzeichens. Also jeweils in der Sprache, wo ich es einfacher fand.

Hier die Unsigned Multiplikation, ist sicherlich noch optimierungsfähig:
Delphi-Quellcode:
procedure MulPInt128U256U(const Result: PInt128; Left, Right: PInt128);
//Result: 256 bit --> zwei PInt128 (Lo, Hi) nacheinander!!!
//Unsigned!
//for 32-bit; nur wenn ss = ds, was bei Delphi immer zu sein scheint
//(kleine Optimierung durch Voraussetzung)
var
  buf1, buf2, buf3: array[0..4] of LongWord;
asm
//mul performs an unsigned multiplication of the operand and the accumulator.
//[...] If the operand is a double word, the processor multiplies it
//by the contents of EAX and returns the 64-bit result in EDX and EAX. mul
//sets CF and OF when the upper half of the result is nonzero, otherwise they
//are cleared. Rules for the operand are the same as for the inc instruction.

  push ebx
  push esi

  push eax
  push ecx

  mov esi, edx

  mov ecx, [ecx]
  mov ebx, eax
  call @@multiply

  mov ebx, [esp] //=original ecx
  mov ecx, [ebx+4]
  lea ebx, [buf1]
  call @@multiply

  mov ebx, [esp]
  mov ecx, [ebx+8]
  lea ebx, [buf2]
  call @@multiply

  mov ebx, [esp]
  mov ecx, [ebx+12]
  lea ebx, [buf3]
  call @@multiply

  //pop ecx
  add esp, 4 //Wert von ecx nicht wichtig, add schneller
  pop eax

  lea ebx, [buf1]
  xor ecx, ecx
  call @@add

  lea ebx, [buf2]
  call @@add

  lea ebx, [buf3]
  call @@add

  jmp @@end

  @@multiply:
    //eine Zeile multiplizieren
    //ecx = dword-faktor
    //abx = buf für ergebnis
    //und esi: Zeigt auf 128-Bit-Faktor!!!

    mov eax, [esi]
    mul ecx
    mov [ebx], eax
    mov [ebx+4], edx

    mov eax, [esi+4]
    mul ecx
    add[ebx+4], eax
    adc edx, 0
    mov [ebx+8], edx
    jc @@m1 //Übertrag speichern
      mov [ebx+12], 0
      jmp @@m2
    @@m1:
      mov [ebx+12], 1
    @@m2:

    mov eax, [esi+8]
    mul ecx
    add [ebx+8], eax
    adc [ebx+12], edx
    jc @@m3 //Übertrag speichern
      mov [ebx+16], 0
      jmp @@m4
    @@m3:
      mov [ebx+16], 1
    @@m4:

    mov eax, [esi+12]
    mul ecx
    add [ebx+12], eax
    adc [ebx+16], edx

  //Wenn jetzt noch ein Bit Übertrag ist, scheiss ich drauf!
  //Oder quäle mich wenn es immer total falsch rechnet :(
  //Meiner Berechnung sollten die 160 Bit (5*32 Bit) aber genau passen für eine
  //Multiplikation 32-Bit mit 128-Bit
  ret

  @@addcarry:
    //eax --> Ergebnis-Stelle
    //ebx ist wieder buf --> Stelle
    //ecx - letzter Übertrag

    add ecx, [ebx] //Übertrag verrechnen

    jnc @@a1 //neuer Übertrag?
      mov edx, 1 //Speichern
    jmp @@a2
    @@a1: //kein neuer Übertrag
      xor edx, edx //Speichern
    @@a2:

    add [eax], ecx //eigentliche Addition

    jnc @@a3 //neuer Übertrag?
      inc edx //Speichern/erhöhen
    @@a3:

    mov ecx, edx
  ret

  @@add:
    //EAX -> Ergebnis, letzte Stelle
    //(nicht sie Stelle, die zuerst addiert werden soll, sondern DAVOR)
    //zwei der Zeilen zusammenaddieren
    //ebx ist wieder buf
    //ecx - letzter Übertrag

    add eax, 4 //ne Stelle weiterrutschen

    push eax //4 Stellen addieren...
      call @@addcarry

      add eax, 4
      add ebx, 4
      call @@addcarry

      add eax, 4
      add ebx, 4
      call @@addcarry

      add eax, 4
      add ebx, 4
      call @@addcarry
    pop eax
  ret

  @@end:

  pop esi
  pop ebx
end;
Und hier die Vorzeichenbehaftete Multiplikation, die darauf zurückgreift:
Delphi-Quellcode:
procedure MulPInt128(const Result: PInt128; Left, Right: PInt128);
var
  Negative: Boolean;
  Help1, Help2: PInt128;
  UseH1, UseH2: Boolean;
  Reslt{$IFOPT Q+}, Reslt2{$ENDIF}: PInt128;
begin
  Help1 := nil; UseH1 := false;
  Help2 := nil; UseH2 := false;

  Negative := false;

  try

    if Left^ < 0 then
    begin
      Negative := true;
      UseH1 := true;
      New(Help1);
      Help1^ := -Left^;
    end
    else
      Help1 := Left;
    if Right^ < 0 then
    begin
      Negative := not Negative;
      UseH2 := true;
      New(Help2);
      Help2^ := -Right^;
    end
    else
      Help2 := Right;

    GetMem(Reslt, SizeOf(Int128) * 2);
    try
      MulPInt128U256U(Reslt, Help1, Help2);

      Result^ := Reslt^;

      {$IFOPT Q+}
      Reslt2 := Reslt;
      inc(Reslt2);

      if (Reslt2^ <> 0) or (Reslt < 0) then
        raise EIntOverflow.Create('Integer Overflow.');
      {$ENDIF}
    finally
      FreeMem(Reslt);
    end;

    if Negative then
      Result^ := -Result^;
  finally
    if UseH2 then Dispose(Help2);
    if UseH1 then Dispose(Help1);
  end;
end;
Mfg
FAlter
Felix Alter
  Mit Zitat antworten Zitat