gcc/testsuite/ada/acats/tests/cxg/cxg2020.a

   1 -- CXG2020.A
   2 --
   3 --                             Grant of Unlimited Rights
   4 --
   5 --     Under contracts F33600-87-D-0337, F33600-84-D-0280, MDA903-79-C-0687,
   6 --     F08630-91-C-0015, and DCA100-97-D-0025, the U.S. Government obtained
   7 --     unlimited rights in the software and documentation contained herein.
   8 --     Unlimited rights are defined in DFAR 252.227-7013(a)(19).  By making
   9 --     this public release, the Government intends to confer upon all
  10 --     recipients unlimited rights  equal to those held by the Government.
  11 --     These rights include rights to use, duplicate, release or disclose the
  12 --     released technical data and computer software in whole or in part, in
  13 --     any manner and for any purpose whatsoever, and to have or permit others
  14 --     to do so.
  15 --
  16 --                                    DISCLAIMER
  17 --
  18 --     ALL MATERIALS OR INFORMATION HEREIN RELEASED, MADE AVAILABLE OR
  19 --     DISCLOSED ARE AS IS.  THE GOVERNMENT MAKES NO EXPRESS OR IMPLIED
  20 --     WARRANTY AS TO ANY MATTER WHATSOEVER, INCLUDING THE CONDITIONS OF THE
  21 --     SOFTWARE, DOCUMENTATION OR OTHER INFORMATION RELEASED, MADE AVAILABLE
  22 --     OR DISCLOSED, OR THE OWNERSHIP, MERCHANTABILITY, OR FITNESS FOR A
  23 --     PARTICULAR PURPOSE OF SAID MATERIAL.
  24 --*
  25 --
  26 -- OBJECTIVE:
  27 --      Check that the complex SQRT function returns
  28 --      a result that is within the error bound allowed.
  29 --
  30 -- TEST DESCRIPTION:
  31 --      This test consists of a generic package that is
  32 --      instantiated to check complex numbers based upon
  33 --      both Float and a long float type.
  34 --      The test for each floating point type is divided into
  35 --      several parts:
  36 --         Special value checks where the result is a known constant.
  37 --         Checks that use an identity for determining the result.
  38 --
  39 -- SPECIAL REQUIREMENTS
  40 --      The Strict Mode for the numerical accuracy must be
  41 --      selected.  The method by which this mode is selected
  42 --      is implementation dependent.
  43 --
  44 -- APPLICABILITY CRITERIA:
  45 --      This test applies only to implementations supporting the
  46 --      Numerics Annex.
  47 --      This test only applies to the Strict Mode for numerical
  48 --      accuracy.
  49 --
  50 --
  51 -- CHANGE HISTORY:
  52 --      24 Mar 96   SAIC    Initial release for 2.1
  53 --      17 Aug 96   SAIC    Incorporated reviewer comments.
  54 --      03 Jun 98   EDS     Added parens to ensure that the expression is not
  55 --                          evaluated by multiplying its two large terms
  56 --                          together and overflowing.
  57 --!
  58
  59 --
  60 -- References:
  61 --
  62 -- W. J. Cody
  63 -- CELEFUNT: A Portable Test Package for Complex Elementary Functions
  64 -- Algorithm 714, Collected Algorithms from ACM.
  65 -- Published in Transactions On Mathematical Software,
  66 -- Vol. 19, No. 1, March, 1993, pp. 1-21.
  67 --
  68 -- CRC Standard Mathematical Tables
  69 -- 23rd Edition
  70 --
  71
  72 with System;
  73 with Report;
  74 with Ada.Numerics.Generic_Complex_Types;
  75 with Ada.Numerics.Generic_Complex_Elementary_Functions;
  76 procedure CXG2020 is
  77    Verbose : constant Boolean := False;
  78    -- Note that Max_Samples is the number of samples taken in
  79    -- both the real and imaginary directions.  Thus, for Max_Samples
  80    -- of 100 the number of values checked is 10000.
  81    Max_Samples : constant := 100;
  82
  83    E  : constant := Ada.Numerics.E;
  84    Pi : constant := Ada.Numerics.Pi;
  85
  86    -- CRC Standard Mathematical Tables;  23rd Edition; pg 738
  87    Sqrt2 : constant :=
  88         1.41421_35623_73095_04880_16887_24209_69807_85696_71875_37695;
  89    Sqrt3 : constant :=
  90         1.73205_08075_68877_29352_74463_41505_87236_69428_05253_81039;
  91
  92    generic
  93       type Real is digits <>;
  94    package Generic_Check is
  95       procedure Do_Test;
  96    end Generic_Check;
  97
  98    package body Generic_Check is
  99       package Complex_Type is new
 100            Ada.Numerics.Generic_Complex_Types (Real);
 101       use Complex_Type;
 102
 103       package CEF is new
 104            Ada.Numerics.Generic_Complex_Elementary_Functions (Complex_Type);
 105
 106       function Sqrt (X : Complex) return Complex renames CEF.Sqrt;
 107
 108       -- flag used to terminate some tests early
 109       Accuracy_Error_Reported : Boolean := False;
 110
 111
 112       procedure Check (Actual, Expected : Real;
 113                        Test_Name : String;
 114                        MRE : Real) is
 115          Max_Error : Real;
 116          Rel_Error : Real;
 117          Abs_Error : Real;
 118       begin
 119          -- In the case where the expected result is very small or 0
 120          -- we compute the maximum error as a multiple of Model_Epsilon
 121          -- instead of Model_Epsilon and Expected.
 122          Rel_Error := MRE * (abs Expected * Real'Model_Epsilon);
 123          Abs_Error := MRE * Real'Model_Epsilon;
 124          if Rel_Error > Abs_Error then
 125             Max_Error := Rel_Error;
 126          else
 127             Max_Error := Abs_Error;
 128          end if;
 129
 130          if abs (Actual - Expected) > Max_Error then
 131             Accuracy_Error_Reported := True;
 132             Report.Failed (Test_Name &
 133                            " actual: " & Real'Image (Actual) &
 134                            " expected: " & Real'Image (Expected) &
 135                            " difference: " & Real'Image (Actual - Expected) &
 136                            " max err:" & Real'Image (Max_Error) );
 137          elsif Verbose then
 138             if Actual = Expected then
 139                Report.Comment (Test_Name & "  exact result");
 140             else
 141                Report.Comment (Test_Name & "  passed");
 142             end if;
 143          end if;
 144       end Check;
 145
 146
 147       procedure Check (Actual, Expected : Complex;
 148                        Test_Name : String;
 149                        MRE : Real) is
 150       begin
 151          Check (Actual.Re, Expected.Re, Test_Name & " real part", MRE);
 152          Check (Actual.Im, Expected.Im, Test_Name & " imaginary part", MRE);
 153       end Check;
 154
 155
 156       procedure Special_Value_Test is
 157          -- In the following tests the expected result is accurate
 158          -- to the machine precision so the minimum guaranteed error
 159          -- bound can be used if the argument is exact.
 160          --
 161          -- One or i is added to the actual and expected results in
 162          -- order to prevent the expected result from having a
 163          -- real or imaginary part of 0.  This is to allow a reasonable
 164          -- relative error for that component.
 165          Minimum_Error : constant := 6.0;
 166          Z1, Z2 : Complex;
 167       begin
 168          Check (Sqrt(9.0+0.0*i) + i,
 169                 3.0+1.0*i,
 170                 "sqrt(9+0i)+i",
 171                 Minimum_Error);
 172          Check (Sqrt (-2.0 + 0.0 * i) + 1.0,
 173                 1.0 + Sqrt2 * i,
 174                 "sqrt(-2)+1 ",
 175                 Minimum_Error);
 176
 177          -- make sure no exception occurs when taking the sqrt of
 178          -- very large and very small values.
 179
 180          Z1 := (Real'Safe_Last * 0.9, Real'Safe_Last * 0.9);
 181          Z2 := Sqrt (Z1);
 182          begin
 183             Check (Z2 * Z2,
 184                    Z1,
 185                    "sqrt((big,big))",
 186                    Minimum_Error + 5.0);  -- +5 for multiply
 187          exception
 188             when others =>
 189                 Report.Failed ("unexpected exception in sqrt((big,big))");
 190          end;
 191
 192          Z1 := (Real'Model_Epsilon * 10.0, Real'Model_Epsilon * 10.0);
 193          Z2 := Sqrt (Z1);
 194          begin
 195             Check (Z2 * Z2,
 196                    Z1,
 197                    "sqrt((little,little))",
 198                    Minimum_Error + 5.0);  -- +5 for multiply
 199          exception
 200             when others =>
 201                 Report.Failed ("unexpected exception in " &
 202                     "sqrt((little,little))");
 203          end;
 204
 205       exception
 206          when Constraint_Error =>
 207             Report.Failed ("Constraint_Error raised in special value test");
 208          when others =>
 209             Report.Failed ("exception in special value test");
 210       end Special_Value_Test;
 211
 212
 213
 214       procedure Exact_Result_Test is
 215          No_Error : constant := 0.0;
 216       begin
 217          -- G.1.2(36);6.0
 218          Check (Sqrt(0.0 + 0.0*i),  0.0 + 0.0 * i, "sqrt(0+0i)", No_Error);
 219
 220          -- G.1.2(37);6.0
 221          Check (Sqrt(1.0 + 0.0*i),  1.0 + 0.0 * i, "sqrt(1+0i)", No_Error);
 222
 223          -- G.1.2(38-39);6.0
 224          Check (Sqrt(-1.0 + 0.0*i),  0.0 + 1.0 * i, "sqrt(-1+0i)", No_Error);
 225
 226          -- G.1.2(40);6.0
 227          if Real'Signed_Zeros then
 228             Check (Sqrt(-1.0-0.0*i), 0.0 - 1.0 * i, "sqrt(-1-0i)", No_Error);
 229          end if;
 230       exception
 231          when Constraint_Error =>
 232             Report.Failed ("Constraint_Error raised in Exact_Result Test");
 233          when others =>
 234             Report.Failed ("exception in Exact_Result Test");
 235       end Exact_Result_Test;
 236
 237
 238       procedure Identity_Test (RA, RB, IA, IB : Real) is
 239       -- Tests an identity over a range of values specified
 240       -- by the 4 parameters.  RA and RB denote the range for the
 241       -- real part while IA and IB denote the range for the
 242       -- imaginary part of the result.
 243       --
 244       -- For this test we use the identity
 245       --    Sqrt(Z*Z) = Z
 246       --
 247
 248          Scale : Real := Real (Real'Machine_Radix) ** (Real'Mantissa / 2 + 4);
 249          W, X, Y, Z : Real;
 250          CX : Complex;
 251          Actual, Expected : Complex;
 252       begin
 253          Accuracy_Error_Reported := False;  -- reset
 254          for II in 1..Max_Samples loop
 255             X :=  (RB - RA) * Real (II) / Real (Max_Samples) + RA;
 256             for J in 1..Max_Samples loop
 257                Y :=  (IB - IA) * Real (J) / Real (Max_Samples) + IA;
 258
 259                -- purify the arguments to minimize roundoff error.
 260                -- We construct the values so that the products X*X,
 261                -- Y*Y, and X*Y are all exact machine numbers.
 262                -- See Cody page 7 and CELEFUNT code.
 263                Z := X * Scale;
 264                W := Z + X;
 265                X := W - Z;
 266                Z := Y * Scale;
 267                W := Z + Y;
 268                Y := W - Z;
 269                  -- G.1.2(21);6.0 - real part of result is non-negative
 270                Expected := Compose_From_Cartesian( abs X,Y);
 271                Z := X*X - Y*Y;
 272                W := X*Y;
 273                CX := Compose_From_Cartesian(Z,W+W);
 274
 275                -- The arguments are now ready so on with the
 276                -- identity computation.
 277                Actual := Sqrt(CX);
 278
 279                Check (Actual, Expected,
 280                       "Identity_1_Test " & Integer'Image (II) &
 281                          Integer'Image (J) & ": Sqrt((" &
 282                          Real'Image (CX.Re) & ", " &
 283                          Real'Image (CX.Im) & ")) ",
 284                       8.5);   -- 6.0 from sqrt, 2.5 from argument.
 285                -- See Cody pg 7-8 for analysis of additional error amount.
 286
 287                if Accuracy_Error_Reported then
 288                  -- only report the first error in this test in order to keep
 289                  -- lots of failures from producing a huge error log
 290                  return;
 291                end if;
 292             end loop;
 293          end loop;
 294
 295       exception
 296          when Constraint_Error =>
 297             Report.Failed
 298                ("Constraint_Error raised in Identity_Test" &
 299                 " for X=(" & Real'Image (X) &
 300                 ", " & Real'Image (X) & ")");
 301          when others =>
 302             Report.Failed ("exception in Identity_Test" &
 303                 " for X=(" & Real'Image (X) &
 304                 ", " & Real'Image (X) & ")");
 305       end Identity_Test;
 306
 307
 308       procedure Do_Test is
 309       begin
 310          Special_Value_Test;
 311          Exact_Result_Test;
 312          -- ranges where the sign is the same and where it
 313          -- differs.
 314          Identity_Test (   0.0,   10.0,       0.0,    10.0);
 315          Identity_Test (   0.0,  100.0,    -100.0,     0.0);
 316       end Do_Test;
 317    end Generic_Check;
 318
 319    -----------------------------------------------------------------------
 320    -----------------------------------------------------------------------
 321    package Float_Check is new Generic_Check (Float);
 322
 323    -- check the floating point type with the most digits
 324    type A_Long_Float is digits System.Max_Digits;
 325    package A_Long_Float_Check is new Generic_Check (A_Long_Float);
 326
 327    -----------------------------------------------------------------------
 328    -----------------------------------------------------------------------
 329
 330
 331 begin
 332    Report.Test ("CXG2020",
 333                 "Check the accuracy of the complex SQRT function");
 334
 335    if Verbose then
 336       Report.Comment ("checking Standard.Float");
 337    end if;
 338
 339    Float_Check.Do_Test;
 340
 341    if Verbose then
 342       Report.Comment ("checking a digits" &
 343                       Integer'Image (System.Max_Digits) &
 344                       " floating point type");
 345    end if;
 346
 347    A_Long_Float_Check.Do_Test;
 348
 349
 350    Report.Result;
 351 end CXG2020;