unit2.pas

用于开发税务票据管理的软件

unit Unit2;

interface
uses
  unit1;
Procedure TRED2(var A:matrx2; N:integer;
                       var D:array of real;var E:array of real);

implementation
Procedure TRED2(var A:matrx2; N:integer;
                       var D:array of real;var E:array of real);
var
     I,J,K,L:integer; H,F,G,SCALE1,ZZ,HH:real;
begin
     If N > 1 Then
     begin
        For I:=N DownTo 2 do
        begin
            L:=I - 1;
            H:=0;
            SCALE1:=0;
            If L > 1 Then 
            Begin 
                For K:=1 To L do
                    SCALE1:=SCALE1 + Abs(A[I, K]);
                If SCALE1 = 0  Then
                    E[I]:=A[I, L]
                else
                begin
                    For K:=1 To L do
                    begin
                        A[I, K]:=A[I, K] / SCALE1;
                        H:=H + Sqr(A[I, K]);
                    end;
                    F:=A[I, L];
                    If F >= 0 then
                        ZZ:=1
                    Else
                        ZZ:=-1;
                    G:=-Sqrt(H) * ZZ;
                    E[I]:=SCALE1 * G;
                    H:=H - F * G;
                    A[I, L]:=F - G;
                    F:=0; 
                    For J:=1 To L do
                    begin
                        A[J, I]:=A[I, J] / H;
                        G:=0;
                        For K:=1 To J do
                            G:=G + A[J, K] * A[I, K];
                        If L > J Then
                        begin
                            For K:=J + 1 To L do
                                G:=G + A[K, J] * A[I, K];
                        End;
                        E[J]:=G / H;
                        F:=F + E[J] * A[I, J];
                    end; 
                    HH:=F / (H + H);
                    For J:=1 To L do
                    begin
                        F:=A[I, J];
                        G:=E[J] - HH * F;
                        E[J]:=G;
                        For K:=1 To J do
                            A[J, K]:=A[J, K] - F * E[K] - G * A[I, K];
                    end; 
                end;
            end
            else
                E[I]:=A[I, L];
            D[I]:=H;
        end;
    end; 
    //Omit following line if finding only eigenvalues.
    D[1]:=0;
    E[1]:=0;
    For I:=1 To N do
    begin
    //Delete lines from here ...
        L:=I - 1;
        If D[I] <> 0  Then
        begin
            For J:=1 To L do
            begin
                G:=0;
                For K:=1 To L do
                    G:=G + A[I, K] * A[K, J];
                For K:=1 To L do
                    A[K, J]:=A[K, J] - G * A[K, I];
            end; 
        end; 
     //... to here when finding only eibenvalues.
        D[I]:=A[I, I];
    //Also delete lines from here ...
        A[I, I]:=1; 
        If L >= 1 Then
        begin
            For J:=1 To L do
            begin
                A[I, J]:=0;
                A[J, I]:=0;
            end;
        end;
    //... to here when finding only eigenvalues.
    end;
end;
end.