axiom-source 20170501-3 / usr / share / axiom-20170501 / src / algebra / D01AGNT.spad

)abbrev package D01AGNT d01AgentsPackage
++ Author: Brian Dupee
++ Date Created: March 1994
++ Date Last Updated: December 1997
++ Description:
++ \axiomType{d01AgentsPackage} is a package of numerical agents to be used
++ to investigate attributes of an input function so as to decide the
++ \axiomFun{measure} of an appropriate numerical integration routine.
++ It contains functions \axiomFun{rangeIsFinite} to test the input range and
++ \axiomFun{functionIsContinuousAtEndPoints} to check for continuity at 
++ the end points of the range.

d01AgentsPackage() : SIG == CODE where

  EF2   ==> ExpressionFunctions2
  EFI  ==> Expression Fraction Integer
  FI  ==> Fraction Integer
  LEDF  ==> List Expression DoubleFloat
  KEDF  ==> Kernel Expression DoubleFloat
  EEDF  ==> Equation Expression DoubleFloat
  EDF  ==> Expression DoubleFloat
  PDF  ==> Polynomial DoubleFloat
  LDF  ==> List DoubleFloat
  SDF  ==> Stream DoubleFloat
  DF  ==> DoubleFloat
  F  ==> Float
  ST  ==> String
  LST  ==> List String
  SI  ==> SingleInteger
  SOCDF  ==> Segment OrderedCompletion DoubleFloat
  OCDF  ==> OrderedCompletion DoubleFloat
  OCEDF  ==> OrderedCompletion Expression DoubleFloat
  EOCEFI  ==> Equation OrderedCompletion Expression Fraction Integer
  OCEFI   ==> OrderedCompletion Expression Fraction Integer
  OCFI    ==> OrderedCompletion Fraction Integer
  NIA  ==> Record(var:Symbol,fn:EDF,range:SOCDF,abserr:DF,relerr:DF)
  INT  ==> Integer
  CTYPE  ==> Union(continuous: "Continuous at the end points",
             lowerSingular: "There is a singularity at the lower end point",
              upperSingular: "There is a singularity at the upper end point",
               bothSingular: "There are singularities at both end points",
                notEvaluated: "End point continuity not yet evaluated")
  RTYPE  ==> Union(finite: "The range is finite",
              lowerInfinite: "The bottom of range is infinite",
                upperInfinite: "The top of range is infinite",
                  bothInfinite: "Both top and bottom points are infinite",
                    notEvaluated: "Range not yet evaluated")
  STYPE  ==> Union(str:SDF,
                   notEvaluated:"Internal singularities not yet evaluated")
  ATT  ==> Record(endPointContinuity:CTYPE,
                    singularitiesStream:STYPE,range:RTYPE)
  ROA  ==> Record(key:NIA,entry:ATT)

  SIG ==> with
    
    rangeIsFinite : NIA -> RTYPE
      ++ rangeIsFinite(args) tests the endpoints of \spad{args.range} for 
      ++ infinite end points. 

    functionIsContinuousAtEndPoints : NIA -> CTYPE
      ++ functionIsContinuousAtEndPoints(args) uses power series limits
      ++ to check for problems at the end points of the range of \spad{args}.

    getlo : SOCDF -> DF
      ++ getlo(x) gets the \axiomType{DoubleFloat} equivalent of
      ++ the first endpoint of the range \axiom{x}

    gethi : SOCDF -> DF
      ++ gethi(x) gets the \axiomType{DoubleFloat} equivalent of
      ++ the second endpoint of the range \axiom{x}

    functionIsOscillatory : NIA -> F
      ++ functionIsOscillatory(a) tests whether the function \spad{a.fn}
      ++ has many zeros of its derivative.

    problemPoints : (EDF, Symbol, SOCDF) -> List DF
      ++ problemPoints(f,var,range) returns a list of possible problem points
      ++ by looking at the zeros of the denominator of the function if it
      ++ can be retracted to \axiomType{Polynomial DoubleFloat}.

    singularitiesOf : NIA -> SDF
      ++ singularitiesOf(args) returns a list of potential 
      ++ singularities of the function within the given range

    df2st : DF -> String 
      ++ df2st(n) coerces a \axiomType{DoubleFloat} to \axiomType{String}

    ldf2lst : LDF -> LST
      ++ ldf2lst(ln) coerces a List of \axiomType{DoubleFloat} to 
      ++ \axiomType{List String}

    sdf2lst : SDF -> LST
      ++ sdf2lst(ln) coerces a Stream of \axiomType{DoubleFloat} to 
      ++ \axiomType{List String}

    commaSeparate : LST -> ST
      ++ commaSeparate(l) produces a comma separated string from a 
      ++ list of strings.

    changeName : (Symbol,Symbol,Result) -> Result
      ++ changeName(s,t,r) changes the name of item \axiom{s} in \axiom{r}
      ++ to \axiom{t}.

  CODE ==> ExpertSystemContinuityPackage add

    import ExpertSystemToolsPackage
    import ExpertSystemContinuityPackage

    -- local functions
    ocdf2ocefi : OCDF -> OCEFI
    rangeOfArgument : (KEDF, NIA) -> DF
    continuousAtPoint? : (EFI,EOCEFI) -> Boolean
    rand:(SOCDF,INT) -> LDF 
    eval:(EDF,Symbol,LDF) -> LDF
    numberOfSignChanges:LDF -> INT
    rangeIsFiniteFunction:NIA -> RTYPE
    functionIsContinuousAtEndPointsFunction:NIA -> CTYPE
 
    changeName(s:Symbol,t:Symbol,r:Result):Result ==
      a := remove!(s,r)$Result
      a case Any =>
        insert!([t,a],r)$Result
        r
      r

    commaSeparate(l:LST):ST ==
      empty?(l)$LST => ""
      (#(l) = 1) => concat(l)$ST
      f := first(l)$LST
      t := [concat([", ",l.i])$ST for i in 2..#(l)]
      concat(f,concat(t)$ST)$ST

    rand(seg:SOCDF,n:INT):LDF ==
      -- produced a sorted list of random numbers in the given range
      l:DF := getlo seg
      s:DF := (gethi seg) - l
      seed:INT := random()$INT
      dseed:DF := seed :: DF
      r:LDF := [(((random(seed)$INT) :: DF)*s/dseed + l) for i in 1..n]
      sort(r)$LDF

    eval(f:EDF,var:Symbol,l:LDF):LDF ==
      empty?(l)$LDF => [0$DF]
      ve := var::EDF
      [retract(eval(f,equation(ve,u::EDF)$EEDF)$EDF)@DF for u in l]

    numberOfSignChanges(l:LDF):INT ==
      -- calculates the number of sign changes in a list
      a := 0$INT
      empty?(l)$LDF => 0
      for i in 2..# l repeat
        if negative?(l.i*l.(i-1))  then
          a := a + 1
      a

    rangeOfArgument(k: KEDF, args:NIA): DF ==
      Args := copy args
      Args.fn := arg := first(argument(k)$KEDF)$LEDF
      functionIsContinuousAtEndPoints(Args) case continuous =>
        r:SOCDF := args.range
        low:EDF := (getlo r) :: EDF
        high:EDF := (gethi r) :: EDF
        eql := equation(a := args.var :: EDF, low)$EEDF
        eqh := equation(a, high)$EEDF
        e1 := (numeric(eval(arg,eql)$EDF)$Numeric(DF)) :: DF
        e2 := (numeric(eval(arg,eqh)$EDF)$Numeric(DF)) :: DF
        e2-e1
      0$DF

    ocdf2ocefi(r:OCDF):OCEFI ==
      finite?(r)$OCDF => (edf2efi(((retract(r)@DF)$OCDF)::EDF))::OCEFI
      r pretend OCEFI

    continuousAtPoint?(f:EFI,e:EOCEFI):Boolean ==
      (l := limit(f,e)$PowerSeriesLimitPackage(FI,EFI)) case OCEFI =>
                       finite?(l :: OCEFI)
      -- if the left hand limit equals the right hand limit, or if neither
      -- side has a limit at this point, the return type of  limit() is
      -- Union(Ordered Completion Expression Fraction Integer,"failed")
      false

    -- exported functions
    
    rangeIsFiniteFunction(args:NIA): RTYPE ==
      -- rangeIsFinite(x) tests the endpoints of x.range for infinite
      -- end points. 
      --             [-inf,  inf]  =>  4
      --             [ x  ,  inf]  =>  3
      --             [-inf,  x  ]  =>  1
      --             [ x  ,  y  ]  =>  0
      fr:SI := (3::SI * whatInfinity(hi(args.range))$OCDF 
                      - whatInfinity(lo(args.range))$OCDF)
      fr = 0 => ["The range is finite"]
      fr = 1 => ["The bottom of range is infinite"]
      fr = 3 => ["The top of range is infinite"]
      fr = 4 => ["Both top and bottom points are infinite"]
      error("rangeIsFinite",["this is not a valid range"])$ErrorFunctions

    rangeIsFinite(args:NIA): RTYPE ==
      nia := copy args
      (t := showAttributes(nia)$IntegrationFunctionsTable) case ATT =>
        s := coerce(t)@ATT
        s.range case notEvaluated => 
          s.range := rangeIsFiniteFunction(nia)
          r:ROA := [nia,s]
          insert!(r)$IntegrationFunctionsTable
          s.range
        s.range
      a:ATT := [["End point continuity not yet evaluated"],
                  ["Internal singularities not yet evaluated"],
                      e:=rangeIsFiniteFunction(nia)]
      r:ROA := [nia,a]
      insert!(r)$IntegrationFunctionsTable
      e

    functionIsContinuousAtEndPointsFunction(args:NIA):CTYPE ==
      v := args.var :: EFI :: OCEFI
      high:OCEFI := ocdf2ocefi(hi(args.range))
      low:OCEFI := ocdf2ocefi(lo(args.range))
      f := edf2efi(args.fn)
      l:Boolean := continuousAtPoint?(f,equation(v,low)$EOCEFI)
      h:Boolean := continuousAtPoint?(f,equation(v,high)$EOCEFI)
      l and h => ["Continuous at the end points"]
      l => ["There is a singularity at the upper end point"]
      h => ["There is a singularity at the lower end point"]
      ["There are singularities at both end points"]

    functionIsContinuousAtEndPoints(args:NIA): CTYPE ==
      nia := copy args
      (t := showAttributes(nia)$IntegrationFunctionsTable) case ATT =>
        s := coerce(t)@ATT
        s.endPointContinuity case notEvaluated => 
          s.endPointContinuity := functionIsContinuousAtEndPointsFunction(nia)
          r:ROA := [nia,s]
          insert!(r)$IntegrationFunctionsTable
          s.endPointContinuity
        s.endPointContinuity
      a:ATT := [e:=functionIsContinuousAtEndPointsFunction(nia),
                 ["Internal singularities not yet evaluated"],
                   ["Range not yet evaluated"]]
      r:ROA := [nia,a]
      insert!(r)$IntegrationFunctionsTable
      e

    functionIsOscillatory(a:NIA):F ==
      args := copy a
      k := tower(numerator args.fn)$EDF
      p:F := pi()$F
      for i in 1..# k repeat
        is?(ker := k.i, sin :: Symbol) => 
          ra := convert(rangeOfArgument(ker,args))@F
          ra > 2*p => return (ra/p)
        is?(ker, cos :: Symbol) => 
          ra := convert(rangeOfArgument(ker,args))@F
          ra > 2*p => return (ra/p)
      l:LDF := rand(args.range,30)
      l := eval(args.fn,args.var,l)
      numberOfSignChanges(l) :: F   

    singularitiesOf(args:NIA):SDF ==
      nia := copy args
      (t := showAttributes(nia)$IntegrationFunctionsTable) case ATT =>
        s:ATT := coerce(t)@ATT
        p:STYPE := s.singularitiesStream
        p case str => p.str
        e:SDF := singularitiesOf(nia.fn,[nia.var],nia.range)
        if not empty?(e) then
          if less?(e,10)$SDF then extend(e,10)$SDF
        s.singularitiesStream := [e]
        r:ROA := [nia,s]
        insert!(r)$IntegrationFunctionsTable
        e
      e:=singularitiesOf(nia.fn,[nia.var],nia.range)
      if not empty?(e) then
        if less?(e,10)$SDF then extend(e,10)$SDF
      a:ATT := [["End point continuity not yet evaluated"],[e],
                          ["Range not yet evaluated"]]
      r:ROA := [nia,a]
      insert!(r)$IntegrationFunctionsTable
      e