diff options
Diffstat (limited to 'gcc-4.8/gcc/ada/a-ngcefu.adb')
| -rw-r--r-- | gcc-4.8/gcc/ada/a-ngcefu.adb | 708 |
1 files changed, 0 insertions, 708 deletions
diff --git a/gcc-4.8/gcc/ada/a-ngcefu.adb b/gcc-4.8/gcc/ada/a-ngcefu.adb deleted file mode 100644 index edcdb5a72..000000000 --- a/gcc-4.8/gcc/ada/a-ngcefu.adb +++ /dev/null @@ -1,708 +0,0 @@ ------------------------------------------------------------------------------- --- -- --- GNAT RUN-TIME COMPONENTS -- --- -- --- ADA.NUMERICS.GENERIC_COMPLEX_ELEMENTARY_FUNCTIONS -- --- -- --- B o d y -- --- -- --- Copyright (C) 1992-2009, Free Software Foundation, Inc. -- --- -- --- GNAT is free software; you can redistribute it and/or modify it under -- --- terms of the GNU General Public License as published by the Free Soft- -- --- ware Foundation; either version 3, or (at your option) any later ver- -- --- sion. GNAT is distributed in the hope that it will be useful, but WITH- -- --- OUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY -- --- or FITNESS FOR A PARTICULAR PURPOSE. -- --- -- --- As a special exception under Section 7 of GPL version 3, you are granted -- --- additional permissions described in the GCC Runtime Library Exception, -- --- version 3.1, as published by the Free Software Foundation. -- --- -- --- You should have received a copy of the GNU General Public License and -- --- a copy of the GCC Runtime Library Exception along with this program; -- --- see the files COPYING3 and COPYING.RUNTIME respectively. If not, see -- --- <http://www.gnu.org/licenses/>. -- --- -- --- GNAT was originally developed by the GNAT team at New York University. -- --- Extensive contributions were provided by Ada Core Technologies Inc. -- --- -- ------------------------------------------------------------------------------- - -with Ada.Numerics.Generic_Elementary_Functions; - -package body Ada.Numerics.Generic_Complex_Elementary_Functions is - - package Elementary_Functions is new - Ada.Numerics.Generic_Elementary_Functions (Real'Base); - use Elementary_Functions; - - PI : constant := 3.14159_26535_89793_23846_26433_83279_50288_41971; - PI_2 : constant := PI / 2.0; - Sqrt_Two : constant := 1.41421_35623_73095_04880_16887_24209_69807_85696; - Log_Two : constant := 0.69314_71805_59945_30941_72321_21458_17656_80755; - - subtype T is Real'Base; - - Epsilon : constant T := 2.0 ** (1 - T'Model_Mantissa); - Square_Root_Epsilon : constant T := Sqrt_Two ** (1 - T'Model_Mantissa); - Inv_Square_Root_Epsilon : constant T := Sqrt_Two ** (T'Model_Mantissa - 1); - Root_Root_Epsilon : constant T := Sqrt_Two ** - ((1 - T'Model_Mantissa) / 2); - Log_Inverse_Epsilon_2 : constant T := T (T'Model_Mantissa - 1) / 2.0; - - Complex_Zero : constant Complex := (0.0, 0.0); - Complex_One : constant Complex := (1.0, 0.0); - Complex_I : constant Complex := (0.0, 1.0); - Half_Pi : constant Complex := (PI_2, 0.0); - - -------- - -- ** -- - -------- - - function "**" (Left : Complex; Right : Complex) return Complex is - begin - if Re (Right) = 0.0 - and then Im (Right) = 0.0 - and then Re (Left) = 0.0 - and then Im (Left) = 0.0 - then - raise Argument_Error; - - elsif Re (Left) = 0.0 - and then Im (Left) = 0.0 - and then Re (Right) < 0.0 - then - raise Constraint_Error; - - elsif Re (Left) = 0.0 and then Im (Left) = 0.0 then - return Left; - - elsif Right = (0.0, 0.0) then - return Complex_One; - - elsif Re (Right) = 0.0 and then Im (Right) = 0.0 then - return 1.0 + Right; - - elsif Re (Right) = 1.0 and then Im (Right) = 0.0 then - return Left; - - else - return Exp (Right * Log (Left)); - end if; - end "**"; - - function "**" (Left : Real'Base; Right : Complex) return Complex is - begin - if Re (Right) = 0.0 and then Im (Right) = 0.0 and then Left = 0.0 then - raise Argument_Error; - - elsif Left = 0.0 and then Re (Right) < 0.0 then - raise Constraint_Error; - - elsif Left = 0.0 then - return Compose_From_Cartesian (Left, 0.0); - - elsif Re (Right) = 0.0 and then Im (Right) = 0.0 then - return Complex_One; - - elsif Re (Right) = 1.0 and then Im (Right) = 0.0 then - return Compose_From_Cartesian (Left, 0.0); - - else - return Exp (Log (Left) * Right); - end if; - end "**"; - - function "**" (Left : Complex; Right : Real'Base) return Complex is - begin - if Right = 0.0 - and then Re (Left) = 0.0 - and then Im (Left) = 0.0 - then - raise Argument_Error; - - elsif Re (Left) = 0.0 - and then Im (Left) = 0.0 - and then Right < 0.0 - then - raise Constraint_Error; - - elsif Re (Left) = 0.0 and then Im (Left) = 0.0 then - return Left; - - elsif Right = 0.0 then - return Complex_One; - - elsif Right = 1.0 then - return Left; - - else - return Exp (Right * Log (Left)); - end if; - end "**"; - - ------------ - -- Arccos -- - ------------ - - function Arccos (X : Complex) return Complex is - Result : Complex; - - begin - if X = Complex_One then - return Complex_Zero; - - elsif abs Re (X) < Square_Root_Epsilon and then - abs Im (X) < Square_Root_Epsilon - then - return Half_Pi - X; - - elsif abs Re (X) > Inv_Square_Root_Epsilon or else - abs Im (X) > Inv_Square_Root_Epsilon - then - return -2.0 * Complex_I * Log (Sqrt ((1.0 + X) / 2.0) + - Complex_I * Sqrt ((1.0 - X) / 2.0)); - end if; - - Result := -Complex_I * Log (X + Complex_I * Sqrt (1.0 - X * X)); - - if Im (X) = 0.0 - and then abs Re (X) <= 1.00 - then - Set_Im (Result, Im (X)); - end if; - - return Result; - end Arccos; - - ------------- - -- Arccosh -- - ------------- - - function Arccosh (X : Complex) return Complex is - Result : Complex; - - begin - if X = Complex_One then - return Complex_Zero; - - elsif abs Re (X) < Square_Root_Epsilon and then - abs Im (X) < Square_Root_Epsilon - then - Result := Compose_From_Cartesian (-Im (X), -PI_2 + Re (X)); - - elsif abs Re (X) > Inv_Square_Root_Epsilon or else - abs Im (X) > Inv_Square_Root_Epsilon - then - Result := Log_Two + Log (X); - - else - Result := 2.0 * Log (Sqrt ((1.0 + X) / 2.0) + - Sqrt ((X - 1.0) / 2.0)); - end if; - - if Re (Result) <= 0.0 then - Result := -Result; - end if; - - return Result; - end Arccosh; - - ------------ - -- Arccot -- - ------------ - - function Arccot (X : Complex) return Complex is - Xt : Complex; - - begin - if abs Re (X) < Square_Root_Epsilon and then - abs Im (X) < Square_Root_Epsilon - then - return Half_Pi - X; - - elsif abs Re (X) > 1.0 / Epsilon or else - abs Im (X) > 1.0 / Epsilon - then - Xt := Complex_One / X; - - if Re (X) < 0.0 then - Set_Re (Xt, PI - Re (Xt)); - return Xt; - else - return Xt; - end if; - end if; - - Xt := Complex_I * Log ((X - Complex_I) / (X + Complex_I)) / 2.0; - - if Re (Xt) < 0.0 then - Xt := PI + Xt; - end if; - - return Xt; - end Arccot; - - -------------- - -- Arccoth -- - -------------- - - function Arccoth (X : Complex) return Complex is - R : Complex; - - begin - if X = (0.0, 0.0) then - return Compose_From_Cartesian (0.0, PI_2); - - elsif abs Re (X) < Square_Root_Epsilon - and then abs Im (X) < Square_Root_Epsilon - then - return PI_2 * Complex_I + X; - - elsif abs Re (X) > 1.0 / Epsilon or else - abs Im (X) > 1.0 / Epsilon - then - if Im (X) > 0.0 then - return (0.0, 0.0); - else - return PI * Complex_I; - end if; - - elsif Im (X) = 0.0 and then Re (X) = 1.0 then - raise Constraint_Error; - - elsif Im (X) = 0.0 and then Re (X) = -1.0 then - raise Constraint_Error; - end if; - - begin - R := Log ((1.0 + X) / (X - 1.0)) / 2.0; - - exception - when Constraint_Error => - R := (Log (1.0 + X) - Log (X - 1.0)) / 2.0; - end; - - if Im (R) < 0.0 then - Set_Im (R, PI + Im (R)); - end if; - - if Re (X) = 0.0 then - Set_Re (R, Re (X)); - end if; - - return R; - end Arccoth; - - ------------ - -- Arcsin -- - ------------ - - function Arcsin (X : Complex) return Complex is - Result : Complex; - - begin - -- For very small argument, sin (x) = x - - if abs Re (X) < Square_Root_Epsilon and then - abs Im (X) < Square_Root_Epsilon - then - return X; - - elsif abs Re (X) > Inv_Square_Root_Epsilon or else - abs Im (X) > Inv_Square_Root_Epsilon - then - Result := -Complex_I * (Log (Complex_I * X) + Log (2.0 * Complex_I)); - - if Im (Result) > PI_2 then - Set_Im (Result, PI - Im (X)); - - elsif Im (Result) < -PI_2 then - Set_Im (Result, -(PI + Im (X))); - end if; - - return Result; - end if; - - Result := -Complex_I * Log (Complex_I * X + Sqrt (1.0 - X * X)); - - if Re (X) = 0.0 then - Set_Re (Result, Re (X)); - - elsif Im (X) = 0.0 - and then abs Re (X) <= 1.00 - then - Set_Im (Result, Im (X)); - end if; - - return Result; - end Arcsin; - - ------------- - -- Arcsinh -- - ------------- - - function Arcsinh (X : Complex) return Complex is - Result : Complex; - - begin - if abs Re (X) < Square_Root_Epsilon and then - abs Im (X) < Square_Root_Epsilon - then - return X; - - elsif abs Re (X) > Inv_Square_Root_Epsilon or else - abs Im (X) > Inv_Square_Root_Epsilon - then - Result := Log_Two + Log (X); -- may have wrong sign - - if (Re (X) < 0.0 and then Re (Result) > 0.0) - or else (Re (X) > 0.0 and then Re (Result) < 0.0) - then - Set_Re (Result, -Re (Result)); - end if; - - return Result; - end if; - - Result := Log (X + Sqrt (1.0 + X * X)); - - if Re (X) = 0.0 then - Set_Re (Result, Re (X)); - elsif Im (X) = 0.0 then - Set_Im (Result, Im (X)); - end if; - - return Result; - end Arcsinh; - - ------------ - -- Arctan -- - ------------ - - function Arctan (X : Complex) return Complex is - begin - if abs Re (X) < Square_Root_Epsilon and then - abs Im (X) < Square_Root_Epsilon - then - return X; - - else - return -Complex_I * (Log (1.0 + Complex_I * X) - - Log (1.0 - Complex_I * X)) / 2.0; - end if; - end Arctan; - - ------------- - -- Arctanh -- - ------------- - - function Arctanh (X : Complex) return Complex is - begin - if abs Re (X) < Square_Root_Epsilon and then - abs Im (X) < Square_Root_Epsilon - then - return X; - else - return (Log (1.0 + X) - Log (1.0 - X)) / 2.0; - end if; - end Arctanh; - - --------- - -- Cos -- - --------- - - function Cos (X : Complex) return Complex is - begin - return - Compose_From_Cartesian - (Cos (Re (X)) * Cosh (Im (X)), - -(Sin (Re (X)) * Sinh (Im (X)))); - end Cos; - - ---------- - -- Cosh -- - ---------- - - function Cosh (X : Complex) return Complex is - begin - return - Compose_From_Cartesian - (Cosh (Re (X)) * Cos (Im (X)), - Sinh (Re (X)) * Sin (Im (X))); - end Cosh; - - --------- - -- Cot -- - --------- - - function Cot (X : Complex) return Complex is - begin - if abs Re (X) < Square_Root_Epsilon and then - abs Im (X) < Square_Root_Epsilon - then - return Complex_One / X; - - elsif Im (X) > Log_Inverse_Epsilon_2 then - return -Complex_I; - - elsif Im (X) < -Log_Inverse_Epsilon_2 then - return Complex_I; - end if; - - return Cos (X) / Sin (X); - end Cot; - - ---------- - -- Coth -- - ---------- - - function Coth (X : Complex) return Complex is - begin - if abs Re (X) < Square_Root_Epsilon and then - abs Im (X) < Square_Root_Epsilon - then - return Complex_One / X; - - elsif Re (X) > Log_Inverse_Epsilon_2 then - return Complex_One; - - elsif Re (X) < -Log_Inverse_Epsilon_2 then - return -Complex_One; - - else - return Cosh (X) / Sinh (X); - end if; - end Coth; - - --------- - -- Exp -- - --------- - - function Exp (X : Complex) return Complex is - EXP_RE_X : constant Real'Base := Exp (Re (X)); - - begin - return Compose_From_Cartesian (EXP_RE_X * Cos (Im (X)), - EXP_RE_X * Sin (Im (X))); - end Exp; - - function Exp (X : Imaginary) return Complex is - ImX : constant Real'Base := Im (X); - - begin - return Compose_From_Cartesian (Cos (ImX), Sin (ImX)); - end Exp; - - --------- - -- Log -- - --------- - - function Log (X : Complex) return Complex is - ReX : Real'Base; - ImX : Real'Base; - Z : Complex; - - begin - if Re (X) = 0.0 and then Im (X) = 0.0 then - raise Constraint_Error; - - elsif abs (1.0 - Re (X)) < Root_Root_Epsilon - and then abs Im (X) < Root_Root_Epsilon - then - Z := X; - Set_Re (Z, Re (Z) - 1.0); - - return (1.0 - (1.0 / 2.0 - - (1.0 / 3.0 - (1.0 / 4.0) * Z) * Z) * Z) * Z; - end if; - - begin - ReX := Log (Modulus (X)); - - exception - when Constraint_Error => - ReX := Log (Modulus (X / 2.0)) - Log_Two; - end; - - ImX := Arctan (Im (X), Re (X)); - - if ImX > PI then - ImX := ImX - 2.0 * PI; - end if; - - return Compose_From_Cartesian (ReX, ImX); - end Log; - - --------- - -- Sin -- - --------- - - function Sin (X : Complex) return Complex is - begin - if abs Re (X) < Square_Root_Epsilon and then - abs Im (X) < Square_Root_Epsilon then - return X; - end if; - - return - Compose_From_Cartesian - (Sin (Re (X)) * Cosh (Im (X)), - Cos (Re (X)) * Sinh (Im (X))); - end Sin; - - ---------- - -- Sinh -- - ---------- - - function Sinh (X : Complex) return Complex is - begin - if abs Re (X) < Square_Root_Epsilon and then - abs Im (X) < Square_Root_Epsilon - then - return X; - - else - return Compose_From_Cartesian (Sinh (Re (X)) * Cos (Im (X)), - Cosh (Re (X)) * Sin (Im (X))); - end if; - end Sinh; - - ---------- - -- Sqrt -- - ---------- - - function Sqrt (X : Complex) return Complex is - ReX : constant Real'Base := Re (X); - ImX : constant Real'Base := Im (X); - XR : constant Real'Base := abs Re (X); - YR : constant Real'Base := abs Im (X); - R : Real'Base; - R_X : Real'Base; - R_Y : Real'Base; - - begin - -- Deal with pure real case, see (RM G.1.2(39)) - - if ImX = 0.0 then - if ReX > 0.0 then - return - Compose_From_Cartesian - (Sqrt (ReX), 0.0); - - elsif ReX = 0.0 then - return X; - - else - return - Compose_From_Cartesian - (0.0, Real'Copy_Sign (Sqrt (-ReX), ImX)); - end if; - - elsif ReX = 0.0 then - R_X := Sqrt (YR / 2.0); - - if ImX > 0.0 then - return Compose_From_Cartesian (R_X, R_X); - else - return Compose_From_Cartesian (R_X, -R_X); - end if; - - else - R := Sqrt (XR ** 2 + YR ** 2); - - -- If the square of the modulus overflows, try rescaling the - -- real and imaginary parts. We cannot depend on an exception - -- being raised on all targets. - - if R > Real'Base'Last then - raise Constraint_Error; - end if; - - -- We are solving the system - - -- XR = R_X ** 2 - Y_R ** 2 (1) - -- YR = 2.0 * R_X * R_Y (2) - -- - -- The symmetric solution involves square roots for both R_X and - -- R_Y, but it is more accurate to use the square root with the - -- larger argument for either R_X or R_Y, and equation (2) for the - -- other. - - if ReX < 0.0 then - R_Y := Sqrt (0.5 * (R - ReX)); - R_X := YR / (2.0 * R_Y); - - else - R_X := Sqrt (0.5 * (R + ReX)); - R_Y := YR / (2.0 * R_X); - end if; - end if; - - if Im (X) < 0.0 then -- halve angle, Sqrt of magnitude - R_Y := -R_Y; - end if; - return Compose_From_Cartesian (R_X, R_Y); - - exception - when Constraint_Error => - - -- Rescale and try again - - R := Modulus (Compose_From_Cartesian (Re (X / 4.0), Im (X / 4.0))); - R_X := 2.0 * Sqrt (0.5 * R + 0.5 * Re (X / 4.0)); - R_Y := 2.0 * Sqrt (0.5 * R - 0.5 * Re (X / 4.0)); - - if Im (X) < 0.0 then -- halve angle, Sqrt of magnitude - R_Y := -R_Y; - end if; - - return Compose_From_Cartesian (R_X, R_Y); - end Sqrt; - - --------- - -- Tan -- - --------- - - function Tan (X : Complex) return Complex is - begin - if abs Re (X) < Square_Root_Epsilon and then - abs Im (X) < Square_Root_Epsilon - then - return X; - - elsif Im (X) > Log_Inverse_Epsilon_2 then - return Complex_I; - - elsif Im (X) < -Log_Inverse_Epsilon_2 then - return -Complex_I; - - else - return Sin (X) / Cos (X); - end if; - end Tan; - - ---------- - -- Tanh -- - ---------- - - function Tanh (X : Complex) return Complex is - begin - if abs Re (X) < Square_Root_Epsilon and then - abs Im (X) < Square_Root_Epsilon - then - return X; - - elsif Re (X) > Log_Inverse_Epsilon_2 then - return Complex_One; - - elsif Re (X) < -Log_Inverse_Epsilon_2 then - return -Complex_One; - - else - return Sinh (X) / Cosh (X); - end if; - end Tanh; - -end Ada.Numerics.Generic_Complex_Elementary_Functions; |