testbdssunit.pas

maths lib with source

unit testbdssunit;
interface
uses Math, Ap, Sysutils, tsort, descriptivestatistics, bdss;

function TestBDSS(Silent : Boolean):Boolean;
function testbdssunit_test_silent():Boolean;
function testbdssunit_test():Boolean;

implementation

procedure Unset2D(var A : TComplex2DArray);forward;
procedure Unset1D(var A : TReal1DArray);forward;
procedure Unset1DI(var A : TInteger1DArray);forward;
procedure TestSortResults(const ASorted : TReal1DArray;
     const P1 : TInteger1DArray;
     const P2 : TInteger1DArray;
     const AOriginal : TReal1DArray;
     N : Integer;
     var WasErrors : Boolean);forward;


(*************************************************************************
Testing BDSS operations
*************************************************************************)
function TestBDSS(Silent : Boolean):Boolean;
var
    N : Integer;
    I : Integer;
    J : Integer;
    Pass : Integer;
    PassCount : Integer;
    MaxN : Integer;
    MaxNQ : Integer;
    A : TReal1DArray;
    A0 : TReal1DArray;
    AT : TReal1DArray;
    P : TReal2DArray;
    Thresholds : TReal1DArray;
    NI : Integer;
    C : TInteger1DArray;
    P1 : TInteger1DArray;
    P2 : TInteger1DArray;
    Ties : TInteger1DArray;
    PT1 : TInteger1DArray;
    PT2 : TInteger1DArray;
    TieCount : Integer;
    C1 : Integer;
    C0 : Integer;
    NC : Integer;
    PAL : Double;
    PBL : Double;
    PAR : Double;
    PBR : Double;
    CVE : Double;
    CVR : Double;
    Info : Integer;
    Threshold : Double;
    TieBuf : TInteger1DArray;
    CntBuf : TInteger1DArray;
    RMS : Double;
    CVRMS : Double;
    WasErrors : Boolean;
    TiesErrors : Boolean;
    Split2Errors : Boolean;
    OptimalSplitKErrors : Boolean;
    SplitKErrors : Boolean;
begin
    WasErrors := False;
    TiesErrors := False;
    Split2Errors := False;
    SplitKErrors := False;
    OptimalSplitKErrors := False;
    MaxN := 100;
    MaxNQ := 49;
    PassCount := 10;
    
    //
    // Test ties
    //
    N:=1;
    while N<=MaxN do
    begin
        Pass:=1;
        while Pass<=PassCount do
        begin
            
            //
            // (probably) untied data, test DSTie
            //
            Unset1DI(P1);
            Unset1DI(P2);
            Unset1DI(PT1);
            Unset1DI(PT2);
            SetLength(A, N-1+1);
            SetLength(A0, N-1+1);
            SetLength(AT, N-1+1);
            I:=0;
            while I<=N-1 do
            begin
                A[I] := 2*RandomReal-1;
                A0[I] := A[I];
                AT[I] := A[I];
                Inc(I);
            end;
            DSTie(A0, N, Ties, TieCount, P1, P2);
            TagSort(AT, N, PT1, PT2);
            I:=0;
            while I<=N-1 do
            begin
                TiesErrors := TiesErrors or (P1[I]<>PT1[I]);
                TiesErrors := TiesErrors or (P2[I]<>PT2[I]);
                Inc(I);
            end;
            TiesErrors := TiesErrors or (TieCount<>N);
            if TieCount=N then
            begin
                I:=0;
                while I<=N do
                begin
                    TiesErrors := TiesErrors or (Ties[I]<>I);
                    Inc(I);
                end;
            end;
            
            //
            // tied data, test DSTie
            //
            Unset1DI(P1);
            Unset1DI(P2);
            Unset1DI(PT1);
            Unset1DI(PT2);
            SetLength(A, N-1+1);
            SetLength(A0, N-1+1);
            SetLength(AT, N-1+1);
            C1 := 0;
            C0 := 0;
            I:=0;
            while I<=N-1 do
            begin
                A[I] := RandomInteger(2);
                if A[I]=0 then
                begin
                    C0 := C0+1;
                end
                else
                begin
                    C1 := C1+1;
                end;
                A0[I] := A[I];
                AT[I] := A[I];
                Inc(I);
            end;
            DSTie(A0, N, Ties, TieCount, P1, P2);
            TagSort(AT, N, PT1, PT2);
            I:=0;
            while I<=N-1 do
            begin
                TiesErrors := TiesErrors or (P1[I]<>PT1[I]);
                TiesErrors := TiesErrors or (P2[I]<>PT2[I]);
                Inc(I);
            end;
            if (C0=0) or (C1=0) then
            begin
                TiesErrors := TiesErrors or (TieCount<>1);
                if TieCount=1 then
                begin
                    TiesErrors := TiesErrors or (Ties[0]<>0);
                    TiesErrors := TiesErrors or (Ties[1]<>N);
                end;
            end
            else
            begin
                TiesErrors := TiesErrors or (TieCount<>2);
                if TieCount=2 then
                begin
                    TiesErrors := TiesErrors or (Ties[0]<>0);
                    TiesErrors := TiesErrors or (Ties[1]<>C0);
                    TiesErrors := TiesErrors or (Ties[2]<>N);
                end;
            end;
            Inc(Pass);
        end;
        Inc(N);
    end;
    
    //
    // split-2
    //
    
    //
    // General tests for different N's
    //
    N:=1;
    while N<=MaxN do
    begin
        SetLength(A, N-1+1);
        SetLength(C, N-1+1);
        
        //
        // one-tie test
        //
        if N mod 2=0 then
        begin
            I:=0;
            while I<=N-1 do
            begin
                A[I] := N;
                C[I] := I mod 2;
                Inc(I);
            end;
            DSOptimalSplit2(A, C, N, Info, Threshold, PAL, PBL, PAR, PBR, CVE);
            if Info<>-3 then
            begin
                Split2Errors := True;
                Inc(N);
                Continue;
            end;
        end;
        
        //
        // two-tie test
        //
        
        //
        // test #1
        //
        if N>1 then
        begin
            I:=0;
            while I<=N-1 do
            begin
                A[I] := I div ((N+1) div 2);
                C[I] := I div ((N+1) div 2);
                Inc(I);
            end;
            DSOptimalSplit2(A, C, N, Info, Threshold, PAL, PBL, PAR, PBR, CVE);
            if Info<>1 then
            begin
                Split2Errors := True;
                Inc(N);
                Continue;
            end;
            Split2Errors := Split2Errors or (AbsReal(Threshold-0.5)>100*MachineEpsilon);
            Split2Errors := Split2Errors or (AbsReal(PAL-1)>100*MachineEpsilon);
            Split2Errors := Split2Errors or (AbsReal(PBL-0)>100*MachineEpsilon);
            Split2Errors := Split2Errors or (AbsReal(PAR-0)>100*MachineEpsilon);
            Split2Errors := Split2Errors or (AbsReal(PBR-1)>100*MachineEpsilon);
        end;
        Inc(N);
    end;
    
    //
    // Special "CREDIT"-test (transparency coefficient)
    //
    N := 110;
    SetLength(A, N-1+1);
    SetLength(C, N-1+1);
    A[0] := 0.000;
    C[0] := 0;
    A[1] := 0.000;
    C[1] := 0;
    A[2] := 0.000;
    C[2] := 0;
    A[3] := 0.000;
    C[3] := 0;
    A[4] := 0.000;
    C[4] := 0;
    A[5] := 0.000;
    C[5] := 0;
    A[6] := 0.000;
    C[6] := 0;
    A[7] := 0.000;
    C[7] := 1;
    A[8] := 0.000;
    C[8] := 0;
    A[9] := 0.000;
    C[9] := 1;
    A[10] := 0.000;
    C[10] := 0;
    A[11] := 0.000;
    C[11] := 0;
    A[12] := 0.000;
    C[12] := 0;
    A[13] := 0.000;
    C[13] := 0;
    A[14] := 0.000;
    C[14] := 0;
    A[15] := 0.000;
    C[15] := 0;
    A[16] := 0.000;
    C[16] := 0;
    A[17] := 0.000;
    C[17] := 0;
    A[18] := 0.000;
    C[18] := 0;
    A[19] := 0.000;
    C[19] := 0;
    A[20] := 0.000;
    C[20] := 0;
    A[21] := 0.000;
    C[21] := 0;
    A[22] := 0.000;
    C[22] := 1;
    A[23] := 0.000;
    C[23] := 0;
    A[24] := 0.000;
    C[24] := 0;
    A[25] := 0.000;
    C[25] := 0;
    A[26] := 0.000;
    C[26] := 0;
    A[27] := 0.000;
    C[27] := 1;
    A[28] := 0.000;
    C[28] := 0;
    A[29] := 0.000;
    C[29] := 1;
    A[30] := 0.000;
    C[30] := 0;
    A[31] := 0.000;
    C[31] := 1;
    A[32] := 0.000;
    C[32] := 0;
    A[33] := 0.000;
    C[33] := 1;
    A[34] := 0.000;
    C[34] := 0;
    A[35] := 0.030;
    C[35] := 0;
    A[36] := 0.030;
    C[36] := 0;
    A[37] := 0.050;
    C[37] := 0;
    A[38] := 0.070;
    C[38] := 1;
    A[39] := 0.110;
    C[39] := 0;
    A[40] := 0.110;
    C[40] := 1;
    A[41] := 0.120;
    C[41] := 0;
    A[42] := 0.130;
    C[42] := 0;
    A[43] := 0.140;
    C[43] := 0;
    A[44] := 0.140;
    C[44] := 0;
    A[45] := 0.140;
    C[45] := 0;
    A[46] := 0.150;
    C[46] := 0;
    A[47] := 0.150;
    C[47] := 0;
    A[48] := 0.170;
    C[48] := 0;
    A[49] := 0.190;
    C[49] := 1;
    A[50] := 0.200;
    C[50] := 0;
    A[51] := 0.200;
    C[51] := 0;
    A[52] := 0.250;
    C[52] := 0;
    A[53] := 0.250;
    C[53] := 0;
    A[54] := 0.260;
    C[54] := 0;
    A[55] := 0.270;
    C[55] := 0;
    A[56] := 0.280;
    C[56] := 0;
    A[57] := 0.310;
    C[57] := 0;
    A[58] := 0.310;
    C[58] := 0;
    A[59] := 0.330;
    C[59] := 0;
    A[60] := 0.330;
    C[60] := 0;
    A[61] := 0.340;
    C[61] := 0;
    A[62] := 0.340;
    C[62] := 0;
    A[63] := 0.370;
    C[63] := 0;
    A[64] := 0.380;
    C[64] := 1;
    A[65] := 0.380;
    C[65] := 0;
    A[66] := 0.410;
    C[66] := 0;
    A[67] := 0.460;
    C[67] := 0;
    A[68] := 0.520;
    C[68] := 0;
    A[69] := 0.530;
    C[69] := 0;
    A[70] := 0.540;
    C[70] := 0;
    A[71] := 0.560;
    C[71] := 0;
    A[72] := 0.560;
    C[72] := 0;
    A[73] := 0.570;
    C[73] := 0;
    A[74] := 0.600;
    C[74] := 0;
    A[75] := 0.600;
    C[75] := 0;
    A[76] := 0.620;
    C[76] := 0;
    A[77] := 0.650;
    C[77] := 0;
    A[78] := 0.660;
    C[78] := 0;
    A[79] := 0.680;
    C[79] := 0;
    A[80] := 0.700;
    C[80] := 0;
    A[81] := 0.750;
    C[81] := 0;
    A[82] := 0.770;
    C[82] := 0;
    A[83] := 0.770;
    C[83] := 0;
    A[84] := 0.770;
    C[84] := 0;
    A[85] := 0.790;
    C[85] := 0;
    A[86] := 0.810;
    C[86] := 0;
    A[87] := 0.840;
    C[87] := 0;
    A[88] := 0.860;
    C[88] := 0;
    A[89] := 0.870;
    C[89] := 0;
    A[90] := 0.890;
    C[90] := 0;
    A[91] := 0.900;
    C[91] := 1;
    A[92] := 0.900;
    C[92] := 0;
    A[93] := 0.910;
    C[93] := 0;
    A[94] := 0.940;
    C[94] := 0;
    A[95] := 0.950;
    C[95] := 0;
    A[96] := 0.952;
    C[96] := 0;
    A[97] := 0.970;
    C[97] := 0;
    A[98] := 0.970;
    C[98] := 0;
    A[99] := 0.980;
    C[99] := 0;
    A[100] := 1.000;
    C[100] := 0;
    A[101] := 1.000;
    C[101] := 0;
    A[102] := 1.000;
    C[102] := 0;
    A[103] := 1.000;
    C[103] := 0;
    A[104] := 1.000;
    C[104] := 0;
    A[105] := 1.020;
    C[105] := 0;
    A[106] := 1.090;
    C[106] := 0;
    A[107] := 1.130;
    C[107] := 0;
    A[108] := 1.840;
    C[108] := 0;
    A[109] := 2.470;
    C[109] := 0;
    DSOptimalSplit2(A, C, N, Info, Threshold, PAL, PBL, PAR, PBR, CVE);
    if Info<>1 then
    begin
        Split2Errors := True;
    end
    else
    begin
        Split2Errors := Split2Errors or (AbsReal(Threshold-0.195)>100*MachineEpsilon);
        Split2Errors := Split2Errors or (AbsReal(PAL-0.80)>0.02);
        Split2Errors := Split2Errors or (AbsReal(PBL-0.20)>0.02);
        Split2Errors := Split2Errors or (AbsReal(PAR-0.97)>0.02);
        Split2Errors := Split2Errors or (AbsReal(PBR-0.03)>0.02);
    end;
    
    //
    // split-2 fast
    //
    
    //
    // General tests for different N's
    //
    N:=1;
    while N<=MaxN do
    begin
        SetLength(A, N-1+1);
        SetLength(C, N-1+1);
        SetLength(TieBuf, N+1);
        SetLength(CntBuf, 3+1);
        
        //
        // one-tie test
        //
        if N mod 2=0 then
        begin
            I:=0;
            while I<=N-1 do
            begin
                A[I] := N;
                C[I] := I mod 2;
                Inc(I);
            end;
            DSOptimalSplit2Fast(A, C, TieBuf, CntBuf, N, 2, 0.00, Info, Threshold, RMS, CVRMS);
            if Info<>-3 then
            begin
                Split2Errors := True;
                Inc(N);
                Continue;
            end;
        end;
        
        //
        // two-tie test
        //
        
        //
        // test #1
        //
        if N>1 then
        begin
            I:=0;
            while I<=N-1 do
            begin
                A[I] := I div ((N+1) div 2);
                C[I] := I div ((N+1) div 2);
                Inc(I);
            end;
            DSOptimalSplit2Fast(A, C, TieBuf, CntBuf, N, 2, 0.00, Info, Threshold, RMS, CVRMS);
            if Info<>1 then
            begin
                Split2Errors := True;
                Inc(N);
                Continue;
            end;
            Split2Errors := Split2Errors or (AbsReal(Threshold-0.5)>100*MachineEpsilon);
            Split2Errors := Split2Errors or (AbsReal(RMS-0)>100*MachineEpsilon);
            if N=2 then
            begin
                Split2Errors := Split2Errors or (AbsReal(CVRMS-0.5)>100*MachineEpsilon);
            end
            else