separate ( Generic_Elementary_Functions )

function KF_Sin( Y1, Y2 : Common_Float ) return Common_Float is

-- On input, Y1 and Y2 are floating point values in Common_Float;
--   These two variables represent the remainder of the reduced argument
--   X = N * (pi/2) + remainder, where |remainder| <= pi/4.
-- On output, a Common_Float value is returned which represents the 
--   approximation of sin( Y1+Y2 ).

   R1, R2, Q                    : Common_Float;

begin

   R1 := Y1;
   R2 := Y2;

-- The following is the core approximation. We approximate
-- sin(Y1+Y2) by an odd polynomial. The case analysis finds both
-- a suitable floating-point type (less expensive to use than 
-- Common_Float) and an appropriate polynomial approximation
-- that will deliver a result accurate enough with respect to
-- Float_Type'Base'Digits. Note that the upper bounds of the
-- cases below (6, 15, 16, 18, 27, and 33) are attributes
-- of predefined floating types of common systems.

   case Float_Type'Base'Digits is

      when 1..6 =>

         declare
            type Working_Float is digits 6;
            R, S, Poly : Working_Float;
         begin
            R := Working_Float( R1 + R2 );
            S := R * R;
            Poly := S*( 0.16666_66269E-0 -
		    S*( 0.83329_37955E-2 -
		    S*( 0.19729_49430E-3 -
		    S*( 0.17898_67484E-5 ))));
            Q := R1 + ( R2 - (R1+R2)*Common_Float( Poly ) );
         end;

      when 7..15 =>

         declare
            type Working_Float is digits 
               (15+System.Max_Digits - abs(15-System.Max_Digits))/2;
            -- this is min( 15, System.Max_Digits ) 
            R, S, Poly : Working_Float;
         begin
            R := Working_Float( R1 + R2 );
            S := R * R;
            Poly := S*( 0.16666_66666_66666_62965E-00 -
		    S*( 0.83333_33333_33216_57487E-02 -
		    S*( 0.19841_26984_00425_50051E-03 -
		    S*( 0.27557_31863_01054_79460E-05 -
                    S*( 0.25051_96268_35442_39751E-07 -
		    S*( 0.16041_33183_44428_30041E-09 -
		    S*( 0.67810_12820_36054_53944E-12 ))))))); 
            Q := R1 + ( R2 - (R1+R2)*Common_Float( Poly ) );
         end;

      when 16 =>

         declare
            type Working_Float is digits 
               (16+System.Max_Digits - abs(16-System.Max_Digits))/2;
            R, S, Poly : Working_Float; 
         begin
            R := Working_Float( R1 + R2 );
            S := R * R;
            Poly := S*( 0.16666_66666_66666_62965E-00 -
		    S*( 0.83333_33333_33216_57487E-02 -
		    S*( 0.19841_26984_00425_50051E-03 -
		    S*( 0.27557_31863_01054_79460E-05 -
                    S*( 0.25051_96268_35442_39751E-07 -
		    S*( 0.16041_33183_44428_30041E-09 -
		    S*( 0.67810_12820_36054_53944E-12 ))))))); 
            Q := R1 + ( R2 - (R1+R2)*Common_Float( Poly ) );
         end;

      when 17..18 =>

         declare
            type Working_Float is digits
               (18+System.Max_Digits - abs(18-System.Max_Digits))/2;
            R, S, Poly : Working_Float;
         begin
            R := Working_Float( R1 + R2 );
            S := R * R;
            Poly := S*( 0.16666_66666_66666_66663_05265E-00 -
		    S*( 0.83333_33333_33333_18353_27331E-02 -
		    S*( 0.19841_26984_12677_42397_98061E-03 -
		    S*( 0.27557_31922_25915_00034_88464E-05 -
                    S*( 0.25052_10788_96438_73197_48920E-07 -
		    S*( 0.16058_94677_96247_63281_41150E-09 -
		    S*( 0.76372_69506_69707_85652_27420E-12 - 
		    S*( 0.23950_00674_30593_53449_39104E-14 ))))))));
            Q := R1 + ( R2 - (R1+R2)*Common_Float( Poly ) );

         end;

      when 19..27 =>

         declare
            type Working_Float is digits
               (27+System.Max_Digits - abs(27-System.Max_Digits))/2;
            R, S, Poly : Working_Float; 
         begin
            R := Working_Float( R1 + R2 );
            S := R * R;
            Poly := S*( 0.16666_66666_66666_66666_66666_66633_0086E+00 -
	       S*( 0.83333_33333_33333_33333_33330_76943_0210E-02 -
	       S*( 0.19841_26984_12698_41269_83450_93722_6766E-03 -
	       S*( 0.27557_31922_39858_90643_73301_78259_2770E-05 -
               S*( 0.25052_10838_54417_12158_37226_81690_8857E-07 -
	       S*( 0.16059_04383_68185_46636_41083_91388_2129E-09 -
	       S*( 0.76471_63730_91063_93957_96779_08149_5627E-12 - 
	       S*( 0.28114_57081_79321_12323_34834_60010_1450E-14 -
	       S*( 0.82204_32293_08044_34287_39069_99726_9855E-17 - 
	       S*( 0.19438_17832_99374_63087_47164_67848_4712E-19 
                   ))))))))));
            Q := R1 + ( R2 - (R1+R2)*Common_Float( Poly ) );
         end;

      when 28..33 =>

         declare
            type Working_Float is digits
               (33+System.Max_Digits - abs(33-System.Max_Digits))/2;
            R, S, Poly : Working_Float; 
         begin
            R := Working_Float( R1 + R2 );
            S := R * R;
            Poly := S*( 0.16666_66666_66666_66666_66666_66666_6666E+00 -
		S*( 0.83333_33333_33333_33333_33333_33333_0274E-02 -
		S*( 0.19841_26984_12698_41269_84126_98311_2600E-03 -
		S*( 0.27557_31922_39858_90652_55730_24207_5442E-05 -
                S*( 0.25052_10838_54417_18775_03589_32984_9329E-07 -
		S*( 0.16059_04383_68216_14589_51327_14046_5961E-09 -
		S*( 0.76471_63731_81981_25083_98331_27914_0063E-12 - 
		S*( 0.28114_57254_34445_07256_92700_63990_6368E-14 -
		S*( 0.82206_35244_69072_20561_89193_97777_7182E-17 - 
		S*( 0.19572_93878_24940_42631_49760_79163_7901E-19 - 
		S*( 0.38680_05154_37551_58372_98445_80258_9098E-22 -
		S*( 0.63830_99586_87078_27791_97302_21097_2554E-25 
                        ))))))))))));
            Q := R1 + ( R2 - (R1+R2)*Common_Float( Poly ) );
         end;

      when others =>

         raise Program_Error;  -- assumption (1) is violated.

   end case;

-- This completes the core approximation.

   return( Q );


end KF_Sin;