File : taylors.ads


with Reps;
use Reps;

generic

  Size: in Positive;
  Dho: in Natural;
  type Scalar is private;

  with function IsNumeric(S: Scalar) return Boolean is <>;
  with function Scal(K: Integer) return Scalar is <>;
  with function Scal(R: Rep) return Scalar is <>;
  with function Ball(R: Rep) return Scalar is <>;
  with function Ball(S: Scalar) return Scalar is <>;
  with function Err0(S: Scalar) return Boolean is <>;
  with function "<"(S1,S2: Scalar) return Boolean is <>;
  with function IntFloor(S: Scalar) return Integer is <>;
  with function IntCeiling(S: Scalar) return Integer is <>;
  with function Approx(S: Scalar) return Rep is <>;
  with function Inf(S: Scalar) return Rep is <>;
  with function Sup(S: Scalar) return Rep is <>;
  with function SupAbs(S: Scalar) return Rep is <>;
  with function "abs"(S: Scalar) return Scalar is <>;
  with function Min(S1,S2: Scalar) return Scalar is <>;
  with function Max(S1,S2: Scalar) return Scalar is <>;
  with function Cap(R: Radius; S: Scalar) return Scalar is <>;
  with function "-"(S: Scalar) return Scalar is <>;

  with function Center(S: Scalar) return Scalar is <>;
  with procedure ResetCenter(S: in out Scalar) is <>;
  with function ResetCenter(S: Scalar) return Scalar is <>;
  with procedure Split(S: in Scalar; S0,SE: out Scalar) is <>;
  with procedure ErrMult(R: in Rep; S: in out Scalar) is <>;
  with function "+"(S,T: Scalar) return Scalar is <>;
  with function "-"(S,T: Scalar) return Scalar is <>;
  with function "*"(R: Rep; S: Scalar) return Scalar is <>;
  with function "/"(S: Scalar; R: Rep) return Scalar is <>;
  with function Sqr(S: Scalar) return Scalar is <>;
  with function "*"(S,T: Scalar) return Scalar is <>;
  with function Inv(S: Scalar) return Scalar is <>;
  with function "/"(S,T: Scalar) return Scalar is <>;
  with function "**"(S: Scalar; I: Integer) return Scalar is <>;

  with function Cosh(S: Scalar) return Scalar is <>;
  with function Sinh(S: Scalar) return Scalar is <>;
  with function Exp(S: Scalar) return Scalar is <>;
  with function Log(S: Scalar) return Scalar is <>;
  with function Lambda(M,N: Integer) return Scalar is <>;

package Taylors is

  pragma Pure;

  subtype Power is Integer range 0 .. Size;
  subtype Degree is Integer range -1 .. Size;
  type Taylor is private;

  type ScalVec is array(Integer range <>) of Scalar;
  type ScalMat is array(Integer range <>, Integer range <>) of Scalar;
  type TayVec is array(Integer range <>) of Taylor;
  type TayMat is array(Integer range <>, Integer range <>) of Taylor;

  function Deg(P: Taylor) return Degree;                                -- Deg
  function EffDeg(P: Taylor) return Degree;                             -- Effective Deg
  procedure SetDeg(D: in Degree; P: in out Taylor);                     -- Set Deg (try to avoid)
  procedure SetZero(P: out Taylor);                                     -- P := 0
  function IsZero(P: Taylor) return Boolean;                            -- True if P=0
  procedure Cleanup(P: in out Taylor);                                  -- Lower Deg if possible
  procedure Copy(P1: in Taylor; P2: out Taylor);                        -- P2 := P1

  function Component(K: Power; P: Taylor) return Scalar;                -- P[K]
  procedure Component(K: in Power; P1: in Taylor; P2: out Taylor);      -- P2[0] := P1[K]
  procedure SetComponent(D: in Power; S: in Scalar; P: in out Taylor);  -- P[D] := S
  procedure ResetComponent(P: in out Taylor;
                     Kmin: in Natural := 0; Kmax: in Integer := Size);  -- Set P.C[Kmin .. Kmax] := 0
  function ErrComp(P: Taylor) return Scalar;                            -- general error
  procedure SetErrComp(S: in Scalar; P: in out Taylor);                 -- general error := |S|
  procedure ResetErrComp(P: in out Taylor);                             -- general error := 0

private

  type Taylor is
    record
      D: Degree;          --- C(D+1 .. Size) considered 0; function=0 if D<0
      E: Scalar;          --- general error term; ignored if using numeric, or if D<0
      C: ScalVec(Power);  --- C(Dho .. Size) may include higher order errors
    end record;

  ScalZero:  constant Scalar := Scal(Zero);
  Trunc: constant Boolean := IsNumeric(ScalZero);

end Taylors;