Mercurial > hg > CbC > CbC_gcc

comparison gcc/testsuite/ada/acats/tests/cxg/cxg2021.a @ 111:04ced10e8804

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
gcc 7
author kono
date Fri, 27 Oct 2017 22:46:09 +0900
parents
children
comparison
equal deleted inserted replaced
68:561a7518be6b 111:04ced10e8804
1 -- CXG2021.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 SIN and COS functions return
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 -- 27 Mar 96 SAIC Initial release for 2.1
53 -- 22 Aug 96 SAIC No longer skips test for systems with
54 -- more than 20 digits of precision.
55 --
56 --!
57
58 --
59 -- References:
60 --
61 -- W. J. Cody
62 -- CELEFUNT: A Portable Test Package for Complex Elementary Functions
63 -- Algorithm 714, Collected Algorithms from ACM.
64 -- Published in Transactions On Mathematical Software,
65 -- Vol. 19, No. 1, March, 1993, pp. 1-21.
66 --
67 -- CRC Standard Mathematical Tables
68 -- 23rd Edition
69 --
70
71 with System;
72 with Report;
73 with Ada.Numerics.Generic_Complex_Types;
74 with Ada.Numerics.Generic_Complex_Elementary_Functions;
75 procedure CXG2021 is
76 Verbose : constant Boolean := False;
77 -- Note that Max_Samples is the number of samples taken in
78 -- both the real and imaginary directions. Thus, for Max_Samples
79 -- of 100 the number of values checked is 10000.
80 Max_Samples : constant := 100;
81
82 E : constant := Ada.Numerics.E;
83 Pi : constant := Ada.Numerics.Pi;
84
85 generic
86 type Real is digits <>;
87 package Generic_Check is
88 procedure Do_Test;
89 end Generic_Check;
90
91 package body Generic_Check is
92 package Complex_Type is new
93 Ada.Numerics.Generic_Complex_Types (Real);
94 use Complex_Type;
95
96 package CEF is new
97 Ada.Numerics.Generic_Complex_Elementary_Functions (Complex_Type);
98
99 function Sin (X : Complex) return Complex renames CEF.Sin;
100 function Cos (X : Complex) return Complex renames CEF.Cos;
101
102 -- flag used to terminate some tests early
103 Accuracy_Error_Reported : Boolean := False;
104
105 -- The following value is a lower bound on the accuracy
106 -- required. It is normally 0.0 so that the lower bound
107 -- is computed from Model_Epsilon. However, for tests
108 -- where the expected result is only known to a certain
109 -- amount of precision this bound takes on a non-zero
110 -- value to account for that level of precision.
111 Error_Low_Bound : Real := 0.0;
112
113 -- the E_Factor is an additional amount added to the Expected
114 -- value prior to computing the maximum relative error.
115 -- This is needed because the error analysis (Cody pg 17-20)
116 -- requires this additional allowance.
117 procedure Check (Actual, Expected : Real;
118 Test_Name : String;
119 MRE : Real;
120 E_Factor : Real := 0.0) is
121 Max_Error : Real;
122 Rel_Error : Real;
123 Abs_Error : Real;
124 begin
125 -- In the case where the expected result is very small or 0
126 -- we compute the maximum error as a multiple of Model_Epsilon instead
127 -- of Model_Epsilon and Expected.
128 Rel_Error := MRE * Real'Model_Epsilon * (abs Expected + E_Factor);
129 Abs_Error := MRE * Real'Model_Epsilon;
130 if Rel_Error > Abs_Error then
131 Max_Error := Rel_Error;
132 else
133 Max_Error := Abs_Error;
134 end if;
135
136 -- take into account the low bound on the error
137 if Max_Error < Error_Low_Bound then
138 Max_Error := Error_Low_Bound;
139 end if;
140
141 if abs (Actual - Expected) > Max_Error then
142 Accuracy_Error_Reported := True;
143 Report.Failed (Test_Name &
144 " actual: " & Real'Image (Actual) &
145 " expected: " & Real'Image (Expected) &
146 " difference: " & Real'Image (Actual - Expected) &
147 " max err:" & Real'Image (Max_Error) &
148 " efactor:" & Real'Image (E_Factor) );
149 elsif Verbose then
150 if Actual = Expected then
151 Report.Comment (Test_Name & " exact result");
152 else
153 Report.Comment (Test_Name & " passed" &
154 " actual: " & Real'Image (Actual) &
155 " expected: " & Real'Image (Expected) &
156 " difference: " & Real'Image (Actual - Expected) &
157 " max err:" & Real'Image (Max_Error) &
158 " efactor:" & Real'Image (E_Factor) );
159 end if;
160 end if;
161 end Check;
162
163
164 procedure Check (Actual, Expected : Complex;
165 Test_Name : String;
166 MRE : Real;
167 R_Factor, I_Factor : Real := 0.0) is
168 begin
169 Check (Actual.Re, Expected.Re, Test_Name & " real part",
170 MRE, R_Factor);
171 Check (Actual.Im, Expected.Im, Test_Name & " imaginary part",
172 MRE, I_Factor);
173 end Check;
174
175
176 procedure Special_Value_Test is
177 -- In the following tests the expected result is accurate
178 -- to the machine precision so the minimum guaranteed error
179 -- bound can be used if the argument is exact.
180 -- Since the argument involves Pi, we must allow for this
181 -- inexact argument.
182 Minimum_Error : constant := 11.0;
183 begin
184 Check (Sin (Pi/2.0 + 0.0*i),
185 1.0 + 0.0*i,
186 "sin(pi/2+0i)",
187 Minimum_Error + 1.0);
188 Check (Cos (Pi/2.0 + 0.0*i),
189 0.0 + 0.0*i,
190 "cos(pi/2+0i)",
191 Minimum_Error + 1.0);
192 exception
193 when Constraint_Error =>
194 Report.Failed ("Constraint_Error raised in special value test");
195 when others =>
196 Report.Failed ("exception in special value test");
197 end Special_Value_Test;
198
199
200
201 procedure Exact_Result_Test is
202 No_Error : constant := 0.0;
203 begin
204 -- G.1.2(36);6.0
205 Check (Sin(0.0 + 0.0*i), 0.0 + 0.0 * i, "sin(0+0i)", No_Error);
206 Check (Cos(0.0 + 0.0*i), 1.0 + 0.0 * i, "cos(0+0i)", No_Error);
207 exception
208 when Constraint_Error =>
209 Report.Failed ("Constraint_Error raised in Exact_Result Test");
210 when others =>
211 Report.Failed ("exception in Exact_Result Test");
212 end Exact_Result_Test;
213
214
215 procedure Identity_Test (RA, RB, IA, IB : Real) is
216 -- Tests an identity over a range of values specified
217 -- by the 4 parameters. RA and RB denote the range for the
218 -- real part while IA and IB denote the range for the
219 -- imaginary part.
220 --
221 -- For this test we use the identity
222 -- Sin(Z) = Sin(Z-W) * Cos(W) + Cos(Z-W) * Sin(W)
223 -- and
224 -- Cos(Z) = Cos(Z-W) * Cos(W) - Sin(Z-W) * Sin(W)
225 --
226
227 X, Y : Real;
228 Z : Complex;
229 W : constant Complex := Compose_From_Cartesian(0.0625, 0.0625);
230 ZmW : Complex; -- Z - W
231 Sin_ZmW,
232 Cos_ZmW : Complex;
233 Actual1, Actual2 : Complex;
234 R_Factor : Real; -- additional real error factor
235 I_Factor : Real; -- additional imaginary error factor
236 Sin_W : constant Complex := (6.2581348413276935585E-2,
237 6.2418588008436587236E-2);
238 -- numeric stability is enhanced by using Cos(W) - 1.0 instead of
239 -- Cos(W) in the computation.
240 Cos_W_m_1 : constant Complex := (-2.5431314180235545803E-6,
241 -3.9062493377261771826E-3);
242
243
244 begin
245 if Real'Digits > 20 then
246 -- constants used here accurate to 20 digits. Allow 1
247 -- additional digit of error for computation.
248 Error_Low_Bound := 0.00000_00000_00000_0001;
249 Report.Comment ("accuracy checked to 19 digits");
250 end if;
251
252 Accuracy_Error_Reported := False; -- reset
253 for II in 0..Max_Samples loop
254 X := (RB - RA) * Real (II) / Real (Max_Samples) + RA;
255 for J in 0..Max_Samples loop
256 Y := (IB - IA) * Real (J) / Real (Max_Samples) + IA;
257
258 Z := Compose_From_Cartesian(X,Y);
259 ZmW := Z - W;
260 Sin_ZmW := Sin (ZmW);
261 Cos_ZmW := Cos (ZmW);
262
263 -- now for the first identity
264 -- Sin(Z) = Sin(Z-W) * Cos(W) + Cos(Z-W) * Sin(W)
265 -- = Sin(Z-W) * (1+(Cos(W)-1)) + Cos(Z-W) * Sin(W)
266 -- = Sin(Z-W) + Sin(Z-W)*(Cos(W)-1) + Cos(Z-W)*Sin(W)
267
268
269 Actual1 := Sin (Z);
270 Actual2 := Sin_ZmW + (Sin_ZmW * Cos_W_m_1 + Cos_ZmW * Sin_W);
271
272 -- The computation of the additional error factors are taken
273 -- from Cody pages 17-20.
274
275 R_Factor := abs (Re (Sin_ZmW) * Re (1.0 - Cos_W_m_1)) +
276 abs (Im (Sin_ZmW) * Im (1.0 - Cos_W_m_1)) +
277 abs (Re (Cos_ZmW) * Re (Sin_W)) +
278 abs (Re (Cos_ZmW) * Re (1.0 - Cos_W_m_1));
279
280 I_Factor := abs (Re (Sin_ZmW) * Im (1.0 - Cos_W_m_1)) +
281 abs (Im (Sin_ZmW) * Re (1.0 - Cos_W_m_1)) +
282 abs (Re (Cos_ZmW) * Im (Sin_W)) +
283 abs (Im (Cos_ZmW) * Re (1.0 - Cos_W_m_1));
284
285 Check (Actual1, Actual2,
286 "Identity_1_Test " & Integer'Image (II) &
287 Integer'Image (J) & ": Sin((" &
288 Real'Image (Z.Re) & ", " &
289 Real'Image (Z.Im) & ")) ",
290 11.0, R_Factor, I_Factor);
291
292 -- now for the second identity
293 -- Cos(Z) = Cos(Z-W) * Cos(W) - Sin(Z-W) * Sin(W)
294 -- = Cos(Z-W) * (1+(Cos(W)-1) - Sin(Z-W) * Sin(W)
295 Actual1 := Cos (Z);
296 Actual2 := Cos_ZmW + (Cos_ZmW * Cos_W_m_1 - Sin_ZmW * Sin_W);
297
298 -- The computation of the additional error factors are taken
299 -- from Cody pages 17-20.
300
301 R_Factor := abs (Re (Sin_ZmW) * Re (Sin_W)) +
302 abs (Im (Sin_ZmW) * Im (Sin_W)) +
303 abs (Re (Cos_ZmW) * Re (1.0 - Cos_W_m_1)) +
304 abs (Im (Cos_ZmW) * Im (1.0 - Cos_W_m_1));
305
306 I_Factor := abs (Re (Sin_ZmW) * Im (Sin_W)) +
307 abs (Im (Sin_ZmW) * Re (Sin_W)) +
308 abs (Re (Cos_ZmW) * Im (1.0 - Cos_W_m_1)) +
309 abs (Im (Cos_ZmW) * Re (1.0 - Cos_W_m_1));
310
311 Check (Actual1, Actual2,
312 "Identity_2_Test " & Integer'Image (II) &
313 Integer'Image (J) & ": Cos((" &
314 Real'Image (Z.Re) & ", " &
315 Real'Image (Z.Im) & ")) ",
316 11.0, R_Factor, I_Factor);
317
318 if Accuracy_Error_Reported then
319 -- only report the first error in this test in order to keep
320 -- lots of failures from producing a huge error log
321 Error_Low_Bound := 0.0; -- reset
322 return;
323 end if;
324 end loop;
325 end loop;
326
327 Error_Low_Bound := 0.0; -- reset
328 exception
329 when Constraint_Error =>
330 Report.Failed
331 ("Constraint_Error raised in Identity_Test" &
332 " for Z=(" & Real'Image (X) &
333 ", " & Real'Image (Y) & ")");
334 when others =>
335 Report.Failed ("exception in Identity_Test" &
336 " for Z=(" & Real'Image (X) &
337 ", " & Real'Image (Y) & ")");
338 end Identity_Test;
339
340
341 procedure Do_Test is
342 begin
343 Special_Value_Test;
344 Exact_Result_Test;
345 -- test regions where sin and cos have the same sign and
346 -- about the same magnitude. This will minimize subtraction
347 -- errors in the identities.
348 -- See Cody page 17.
349 Identity_Test (0.0625, 10.0, 0.0625, 10.0);
350 Identity_Test ( 16.0, 17.0, 16.0, 17.0);
351 end Do_Test;
352 end Generic_Check;
353
354 -----------------------------------------------------------------------
355 -----------------------------------------------------------------------
356 package Float_Check is new Generic_Check (Float);
357
358 -- check the floating point type with the most digits
359 type A_Long_Float is digits System.Max_Digits;
360 package A_Long_Float_Check is new Generic_Check (A_Long_Float);
361
362 -----------------------------------------------------------------------
363 -----------------------------------------------------------------------
364
365
366 begin
367 Report.Test ("CXG2021",
368 "Check the accuracy of the complex SIN and COS functions");
369
370 if Verbose then
371 Report.Comment ("checking Standard.Float");
372 end if;
373
374 Float_Check.Do_Test;
375
376 if Verbose then
377 Report.Comment ("checking a digits" &
378 Integer'Image (System.Max_Digits) &
379 " floating point type");
380 end if;
381
382 A_Long_Float_Check.Do_Test;
383
384
385 Report.Result;
386 end CXG2021;