--Copyright (c) 1991-2002, The Numerical ALgorithms Group Ltd.
--All rights reserved.
--
--Redistribution and use in source and binary forms, with or without
--modification, are permitted provided that the following conditions are
--met:
--
--    - Redistributions of source code must retain the above copyright
--      notice, this list of conditions and the following disclaimer.
--
--    - Redistributions in binary form must reproduce the above copyright
--      notice, this list of conditions and the following disclaimer in
--      the documentation and/or other materials provided with the
--      distribution.
--
--    - Neither the name of The Numerical ALgorithms Group Ltd. nor the
--      names of its contributors may be used to endorse or promote products
--      derived from this software without specific prior written permission.
--
--THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
--IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
--TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
--PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
--OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
--EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
--PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
--PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
--LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
--NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
--SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

)abbrev category UPXSCCA UnivariatePuiseuxSeriesConstructorCategory
++ Author: Clifton J. Williamson
++ Date Created: 6 February 1990
++ Date Last Updated: 22 March 1990
++ Basic Operations:
++ Related Domains:
++ Also See:
++ AMS Classifications:
++ Keywords: series, Puiseux, Laurent
++ Examples:
++ References:
++ Description:
++   This is a category of univariate Puiseux series constructed
++   from univariate Laurent series.  A Puiseux series is represented
++   by a pair \spad{[r,f(x)]}, where r is a positive rational number and
++   \spad{f(x)} is a Laurent series.  This pair represents the Puiseux
++   series \spad{f(x^r)}.
UnivariatePuiseuxSeriesConstructorCategory(Coef,ULS):_
 Category == Definition where
  Coef : Ring
  ULS  : UnivariateLaurentSeriesCategory Coef
  I  ==> Integer
  RN ==> Fraction Integer

  Definition ==> Join(UnivariatePuiseuxSeriesCategory(Coef),_
                      RetractableTo ULS,CoercibleFrom ULS) with
    puiseux: (RN,ULS) -> %
      ++ \spad{puiseux(r,f(x))} returns \spad{f(x^r)}.
    rationalPower: % -> RN
      ++ \spad{rationalPower(f(x))} returns r where the Puiseux series
      ++ \spad{f(x) = g(x^r)}.
    laurentRep  : % -> ULS
      ++ \spad{laurentRep(f(x))} returns \spad{g(x)} where the Puiseux series
      ++ \spad{f(x) = g(x^r)} is represented by \spad{[r,g(x)]}.
    degree: % -> RN
      ++ \spad{degree(f(x))} returns the degree of the leading term of the
      ++ Puiseux series \spad{f(x)}, which may have zero as a coefficient.
    laurent: % -> ULS
      ++ \spad{laurent(f(x))} converts the Puiseux series \spad{f(x)} to a
      ++ Laurent series if possible. Error: if this is not possible.
    laurentIfCan: % -> Union(ULS,"failed")
      ++ \spad{laurentIfCan(f(x))} converts the Puiseux series \spad{f(x)}
      ++ to a Laurent series if possible.
      ++ If this is not possible, "failed" is returned.

   add

     zero? x == zero? laurentRep x
     retract(x:%):ULS == laurent x
     retractIfCan(x:%):Union(ULS,"failed") == laurentIfCan x

)abbrev domain UPXSCONS UnivariatePuiseuxSeriesConstructor
++ Author: Clifton J. Williamson
++ Date Created: 9 May 1989
++ Date Last Updated: 30 November 1994
++ Basic Operations:
++ Related Domains:
++ Also See:
++ AMS Classifications:
++ Keywords: series, Puiseux, Laurent
++ Examples:
++ References:
++ Description:
++   This package enables one to construct a univariate Puiseux series
++   domain from a univariate Laurent series domain. Univariate
++   Puiseux series are represented by a pair \spad{[r,f(x)]}, where r is
++   a positive rational number and \spad{f(x)} is a Laurent series.
++   This pair represents the Puiseux series \spad{f(x^r)}.

UnivariatePuiseuxSeriesConstructor(Coef,ULS):_
 Exports == Implementation where
  Coef  : Ring
  ULS   : UnivariateLaurentSeriesCategory Coef
  I     ==> Integer
  L     ==> List
  NNI   ==> NonNegativeInteger
  OUT   ==> OutputForm
  PI    ==> PositiveInteger
  RN    ==> Fraction Integer
  ST    ==> Stream Coef
  LTerm ==> Record(k:I,c:Coef)
  PTerm ==> Record(k:RN,c:Coef)
  ST2LP ==> StreamFunctions2(LTerm,PTerm)
  ST2PL ==> StreamFunctions2(PTerm,LTerm)

  Exports ==> UnivariatePuiseuxSeriesConstructorCategory(Coef,ULS)

  Implementation ==> add

--% representation

    Rep := Record(expon:RN,lSeries:ULS)

    getExpon: % -> RN
    getULS  : % -> ULS

    getExpon pxs == pxs.expon
    getULS   pxs == pxs.lSeries

--% creation and destruction

    puiseux(n,ls)   == [n,ls]
    laurentRep x    == getULS x
    rationalPower x == getExpon x
    degree x        == getExpon(x) * degree(getULS(x))

    0 == puiseux(1,0)
    1 == puiseux(1,1)

    monomial(c,k) ==
      k = 0 => c :: %
      negative? k => puiseux(-k,monomial(c,-1))
      puiseux(k,monomial(c,1))

    coerce(ls:ULS) == puiseux(1,ls)
    coerce(r:Coef) == r :: ULS  :: %
    coerce(i:I)    == i :: Coef :: %

    laurentIfCan upxs ==
      r := getExpon upxs
      one? denom r =>
        multiplyExponents(getULS upxs,numer(r) :: PI)
      "failed"

    laurent upxs ==
      (uls := laurentIfCan upxs) case "failed" =>
        error "laurent: Puiseux series has fractional powers"
      uls :: ULS

    multExp: (RN,LTerm) -> PTerm
    multExp(r,lTerm) == [r * lTerm.k,lTerm.c]

    terms upxs == map(multExp(getExpon upxs,#1),terms getULS upxs)$ST2LP

    clearDen: (I,PTerm) -> LTerm
    clearDen(n,lTerm) ==
      (int := retractIfCan(n * lTerm.k)@Union(I,"failed")) case "failed" =>
        error "series: inappropriate denominator"
      [int :: I,lTerm.c]

    series(n,stream) ==
      str := map(clearDen(n,#1),stream)$ST2PL
      puiseux(1/n,series str)

--% normalizations

    rewrite:(%,PI) -> %
    rewrite(upxs,m) ==
      -- rewrites a series in x**r as a series in x**(r/m)
      puiseux((getExpon upxs)*(1/m),multiplyExponents(getULS upxs,m))

    ratGcd: (RN,RN) -> RN
    ratGcd(r1,r2) ==
      -- if r1 = prod(p prime,p ** ep(r1)) and
      -- if r2 = prod(p prime,p ** ep(r2)), then
      -- ratGcd(r1,r2) = prod(p prime,p ** min(ep(r1),ep(r2)))
      gcd(numer r1,numer r2) / lcm(denom r1,denom r2)

    withNewExpon:(%,RN) -> %
    withNewExpon(upxs,r) ==
      rewrite(upxs,numer(getExpon(upxs)/r) pretend PI)

--% predicates

    upxs1 = upxs2 ==
      r1 := getExpon upxs1; r2 := getExpon upxs2
      ls1 := getULS upxs1; ls2 := getULS upxs2
      (r1 = r2) => (ls1 = ls2)
      r := ratGcd(r1,r2)
      m1 := numer(getExpon(upxs1)/r) pretend PI
      m2 := numer(getExpon(upxs2)/r) pretend PI
      multiplyExponents(ls1,m1) = multiplyExponents(ls2,m2)

    pole? upxs == pole? getULS upxs

--% arithmetic

    applyFcn:((ULS,ULS) -> ULS,%,%) -> %
    applyFcn(op,pxs1,pxs2) ==
      r1 := getExpon pxs1; r2 := getExpon pxs2
      ls1 := getULS pxs1; ls2 := getULS pxs2
      (r1 = r2) => puiseux(r1,op(ls1,ls2))
      r := ratGcd(r1,r2)
      m1 := numer(getExpon(pxs1)/r) pretend PI
      m2 := numer(getExpon(pxs2)/r) pretend PI
      puiseux(r,op(multiplyExponents(ls1,m1),multiplyExponents(ls2,m2)))

    pxs1 + pxs2     == applyFcn(#1 +$ULS #2,pxs1,pxs2)
    pxs1 - pxs2     == applyFcn(#1 -$ULS #2,pxs1,pxs2)
    pxs1:% * pxs2:% == applyFcn(#1 *$ULS #2,pxs1,pxs2)

    pxs:% ** n:NNI == puiseux(getExpon pxs,getULS(pxs)**n)

    recip pxs ==
      rec := recip getULS pxs
      rec case "failed" => "failed"
      puiseux(getExpon pxs,rec :: ULS)

    RATALG : Boolean := Coef has Algebra(Fraction Integer)

    elt(upxs1:%,upxs2:%) ==
      uls1 := laurentRep upxs1; uls2 := laurentRep upxs2
      r1 := rationalPower upxs1; r2 := rationalPower upxs2
      (n := retractIfCan(r1)@Union(Integer,"failed")) case Integer =>
        puiseux(r2,uls1(uls2 ** r1))
      RATALG =>
        if zero? (coef := coefficient(uls2,deg := degree uls2)) then
          deg := order(uls2,deg + 1000)
          zero? (coef := coefficient(uls2,deg)) =>
            error "elt: series with many leading zero coefficients"
        -- a fractional power of a Laurent series may not be defined:
        -- if f(x) = c * x**n + ..., then f(x) ** (p/q) will be defined
        -- only if q divides n
        b := lcm(denom r1,deg); c := b quo deg
        mon : ULS := monomial(1,c)
        uls2 := elt(uls2,mon) ** r1
        puiseux(r2*(1/c),elt(uls1,uls2))
      error "elt: rational powers not available for this coefficient domain"

    if Coef has "**": (Coef,Integer) -> Coef and
       Coef has "**": (Coef, RN) -> Coef then
         eval(upxs:%,a:Coef) == eval(getULS upxs,a ** getExpon(upxs))

    if Coef has Field then

      pxs1:% / pxs2:% == applyFcn(#1 /$ULS #2,pxs1,pxs2)

      inv upxs ==
        (invUpxs := recip upxs) case "failed" =>
          error "inv: multiplicative inverse does not exist"
        invUpxs :: %

--% values

    variable upxs == variable getULS upxs
    center   upxs == center   getULS upxs

    coefficient(upxs,rn) ==
      one? denom(n := rn / getExpon upxs) =>
        coefficient(getULS upxs,numer n)
      0

    elt(upxs:%,rn:RN) == coefficient(upxs,rn)

--% other functions

    roundDown: RN -> I
    roundDown rn ==
      -- returns the largest integer <= rn
      (den := denom rn) = 1 => numer rn
      n := (num := numer rn) quo den
      positive?(num) => n
      n - 1

    roundUp: RN -> I
    roundUp rn ==
      -- returns the smallest integer >= rn
      (den := denom rn) = 1 => numer rn
      n := (num := numer rn) quo den
      positive?(num) => n + 1
      n

    order upxs == getExpon upxs * order getULS upxs
    order(upxs,r) ==
      e := getExpon upxs
      ord := order(getULS upxs, n := roundDown(r / e))
      ord = n => r
      ord * e

    truncate(upxs,r) ==
      e := getExpon upxs
      puiseux(e,truncate(getULS upxs,roundDown(r / e)))

    truncate(upxs,r1,r2) ==
      e := getExpon upxs
      puiseux(e,truncate(getULS upxs,roundUp(r1 / e),roundDown(r2 / e)))

    complete upxs == puiseux(getExpon upxs,complete getULS upxs)
    extend(upxs,r) ==
      e := getExpon upxs
      puiseux(e,extend(getULS upxs,roundDown(r / e)))

    map(fcn,upxs) == puiseux(getExpon upxs,map(fcn,getULS upxs))

    characteristic == characteristic$Coef

    -- multiplyCoefficients(f,upxs) ==
      -- r := getExpon upxs
      -- puiseux(r,multiplyCoefficients(f(#1 * r),getULS upxs))

    multiplyExponents(upxs:%,n:RN) ==
      puiseux(n * getExpon(upxs),getULS upxs)
    multiplyExponents(upxs:%,n:PI) ==
      puiseux(n * getExpon(upxs),getULS upxs)

    if Coef has "*": (Fraction Integer, Coef) -> Coef then

      differentiate upxs ==
        r := getExpon upxs
        puiseux(r,differentiate getULS upxs) * monomial(r :: Coef,r-1)

      if Coef has PartialDifferentialRing(Symbol) then

        differentiate(upxs:%,s:Symbol) ==
          (s = variable(upxs)) => differentiate upxs
          dcds := differentiate(center upxs,s)
          map(differentiate(#1,s),upxs) - dcds*differentiate(upxs)

    if Coef has Algebra Fraction Integer then

      coerce(r:RN) == r :: Coef :: %

      ratInv: RN -> Coef
      ratInv r ==
        zero? r => 1
        inv(r) :: Coef

      integrate upxs ==
        not zero? coefficient(upxs,-1) =>
          error "integrate: series has term of order -1"
        r := getExpon upxs
        uls := getULS upxs
        uls := multiplyCoefficients(ratInv(#1 * r + 1),uls)
        monomial(1,1) * puiseux(r,uls)

      if Coef has integrate: (Coef,Symbol) -> Coef and _
         Coef has variables: Coef -> List Symbol then

        integrate(upxs:%,s:Symbol) ==
          (s = variable(upxs)) => integrate upxs
          not entry?(s,variables center upxs) => map(integrate(#1,s),upxs)
          error "integrate: center is a function of variable of integration"

      if Coef has TranscendentalFunctionCategory and _
         Coef has PrimitiveFunctionCategory and _
         Coef has AlgebraicallyClosedFunctionSpace Integer then

        integrateWithOneAnswer: (Coef,Symbol) -> Coef
        integrateWithOneAnswer(f,s) ==
          res := integrate(f,s)$FunctionSpaceIntegration(I,Coef)
          res case Coef => res :: Coef
          first(res :: List Coef)

        integrate(upxs:%,s:Symbol) ==
          (s = variable(upxs)) => integrate upxs
          not entry?(s,variables center upxs) =>
            map(integrateWithOneAnswer(#1,s),upxs)
          error "integrate: center is a function of variable of integration"

      if Coef has Field then
         (upxs:%) ** (q:RN) ==
           num := numer q; den := denom q
           one? den => upxs ** num
           r := rationalPower upxs; uls := laurentRep upxs
           deg := degree uls
           if zero?(coef := coefficient(uls,deg)) then
             deg := order(uls,deg + 1000)
             zero?(coef := coefficient(uls,deg)) =>
               error "power of series with many leading zero coefficients"
           ulsPow := (uls * monomial(1,-deg)$ULS) ** q
           puiseux(r,ulsPow) * monomial(1,deg*q*r)

      applyUnary: (ULS -> ULS,%) -> %
      applyUnary(fcn,upxs) ==
        puiseux(rationalPower upxs,fcn laurentRep upxs)

      exp upxs   == applyUnary(exp,upxs)
      log upxs   == applyUnary(log,upxs)
      sin upxs   == applyUnary(sin,upxs)
      cos upxs   == applyUnary(cos,upxs)
      tan upxs   == applyUnary(tan,upxs)
      cot upxs   == applyUnary(cot,upxs)
      sec upxs   == applyUnary(sec,upxs)
      csc upxs   == applyUnary(csc,upxs)
      asin upxs  == applyUnary(asin,upxs)
      acos upxs  == applyUnary(acos,upxs)
      atan upxs  == applyUnary(atan,upxs)
      acot upxs  == applyUnary(acot,upxs)
      asec upxs  == applyUnary(asec,upxs)
      acsc upxs  == applyUnary(acsc,upxs)
      sinh upxs  == applyUnary(sinh,upxs)
      cosh upxs  == applyUnary(cosh,upxs)
      tanh upxs  == applyUnary(tanh,upxs)
      coth upxs  == applyUnary(coth,upxs)
      sech upxs  == applyUnary(sech,upxs)
      csch upxs  == applyUnary(csch,upxs)
      asinh upxs == applyUnary(asinh,upxs)
      acosh upxs == applyUnary(acosh,upxs)
      atanh upxs == applyUnary(atanh,upxs)
      acoth upxs == applyUnary(acoth,upxs)
      asech upxs == applyUnary(asech,upxs)
      acsch upxs == applyUnary(acsch,upxs)

)abbrev domain UPXS UnivariatePuiseuxSeries
++ Author: Clifton J. Williamson
++ Date Created: 28 January 1990
++ Date Last Updated: June 18, 2010
++ Basic Operations:
++ Related Domains:
++ Also See:
++ AMS Classifications:
++ Keywords: series, Puiseux
++ Examples:
++ References:
++ Description: Dense Puiseux series in one variable
++   \spadtype{UnivariatePuiseuxSeries} is a domain representing Puiseux
++   series in one variable with coefficients in an arbitrary ring.  The
++   parameters of the type specify the coefficient ring, the power series
++   variable, and the center of the power series expansion.  For example,
++   \spad{UnivariatePuiseuxSeries(Integer,x,3)} represents Puiseux series in
++   \spad{(x - 3)} with \spadtype{Integer} coefficients.
UnivariatePuiseuxSeries(Coef,var,cen): Exports == Implementation where
  Coef : Ring
  var  : Symbol
  cen  : Coef
  I   ==> Integer
  L   ==> List
  NNI ==> NonNegativeInteger
  OUT ==> OutputForm
  RN  ==> Fraction Integer
  ST  ==> Stream Coef
  UTS ==> UnivariateTaylorSeries(Coef,var,cen)
  ULS ==> UnivariateLaurentSeries(Coef,var,cen)

  Exports ==> Join(UnivariatePuiseuxSeriesConstructorCategory(Coef,ULS),_
                   PartialDifferentialDomain(%,Variable var),_
                   RetractableTo UTS,CoercibleFrom Variable var) with
    if Coef has Algebra Fraction Integer then
      integrate: (%,Variable(var)) -> %
        ++ \spad{integrate(f(x))} returns an anti-derivative of the power
        ++ series \spad{f(x)} with constant coefficient 0.
        ++ We may integrate a series when we can divide coefficients
        ++ by integers.

  Implementation ==> UnivariatePuiseuxSeriesConstructor(Coef,ULS) add

    Rep := Record(expon:RN,lSeries:ULS)

    getExpon: % -> RN
    getExpon pxs == pxs.expon

    variable upxs == var
    center   upxs == cen

    coerce(uts:UTS) == uts :: ULS :: %

    retractIfCan(upxs:%):Union(UTS,"failed") ==
      (ulsIfCan := retractIfCan(upxs)@Union(ULS,"failed")) case "failed" =>
        "failed"
      retractIfCan(ulsIfCan :: ULS)

    --retract(upxs:%):UTS ==
      --(ulsIfCan := retractIfCan(upxs)@Union(ULS,"failed")) case "failed" =>
        --error "retractIfCan: series has fractional exponents"
      --utsIfCan := retractIfCan(ulsIfCan :: ULS)@Union(UTS,"failed")
      --utsIfCan case "failed" =>
        --error "retractIfCan: series has negative exponents"
      --utsIfCan :: UTS

    coerce(v:Variable(var)) ==
      zero? cen => monomial(1,1)
      monomial(1,1) + monomial(cen,0)

    if Coef has "*": (Fraction Integer, Coef) -> Coef then
      differentiate(upxs:%,v:Variable(var)) == differentiate upxs

    if Coef has Algebra Fraction Integer then
      integrate(upxs:%,v:Variable(var)) == integrate upxs

    if Coef has coerce: Symbol -> Coef then
      if Coef has "**": (Coef,RN) -> Coef then

        roundDown: RN -> I
        roundDown rn ==
          -- returns the largest integer <= rn
          (den := denom rn) = 1 => numer rn
          n := (num := numer rn) quo den
          positive?(num) => n
          n - 1

        stToCoef: (ST,Coef,NNI,NNI) -> Coef
        stToCoef(st,term,n,n0) ==
          (n > n0) or (empty? st) => 0
          frst(st) * term ** n + stToCoef(rst st,term,n + 1,n0)

        approximateLaurent: (ULS,Coef,I) -> Coef
        approximateLaurent(x,term,n) ==
          negative?(m := n - (e := degree x)) => 0
          app := stToCoef(coefficients taylorRep x,term,0,m :: NNI)
          zero? e => app
          app * term ** (e :: RN)

        approximate(x,r) ==
          e := rationalPower(x)
          term := ((variable(x) :: Coef) - center(x)) ** e
          approximateLaurent(laurentRep x,term,roundDown(r / e))

    termOutput:(RN,Coef,OUT) -> OUT
    termOutput(k,c,vv) ==
    -- creates a term c * vv ** k
      k = 0 => c :: OUT
      mon :=
        k = 1 => vv
        vv ** (k :: OUT)
      c = 1 => mon
      c = -1 => -mon
      (c :: OUT) * mon

    showAll?:() -> Boolean
    -- check a global Lisp variable
    showAll?() == true

    termsToOutputForm:(RN,RN,ST,OUT) -> OUT
    termsToOutputForm(m,rat,uu,xxx) ==
      l : L OUT := empty()
      empty? uu => 0@Coef :: OUT
      count : NNI := _$streamCount$Lisp
      n : NNI := 0
      while n <= count and not empty? uu repeat
        if frst(uu) ~= 0 then
          l := concat(termOutput((n :: I) * rat + m,frst uu,xxx),l)
        uu := rst uu
        n := n + 1
      if showAll?() then
        n := count + 1
        while explicitEntries? uu and _
               not eq?(uu,rst uu) repeat
          if frst(uu) ~= 0 then
            l := concat(termOutput((n :: I) * rat + m,frst uu,xxx),l)
          uu := rst uu
          n := n + 1
      l :=
        explicitlyEmpty? uu => l
        eq?(uu,rst uu) and frst uu = 0 => l
        concat(prefix("O" :: OUT,[xxx ** (((n::I) * rat + m) :: OUT)]),l)
      empty? l => 0@Coef :: OUT
      reduce("+",reverse! l)

    coerce(upxs:%):OUT ==
      rat := getExpon upxs; uls := laurentRep upxs
      count : I := _$streamCount$Lisp
      uls := removeZeroes(_$streamCount$Lisp,uls)
      m : RN := (degree uls) * rat
      p := coefficients taylorRep uls
      xxx :=
        zero? cen => var :: OUT
        paren(var :: OUT - cen :: OUT)
      termsToOutputForm(m,rat,p,xxx)

)abbrev package UPXS2 UnivariatePuiseuxSeriesFunctions2
++ Mapping package for univariate Puiseux series
++ Author: Scott C. Morrison
++ Date Created: 5 April 1991
++ Date Last Updated: 5 April 1991
++ Keywords: Puiseux series, map
++ Examples:
++ References:
++ Description:
++   Mapping package for univariate Puiseux series.
++   This package allows one to apply a function to the coefficients of
++   a univariate Puiseux series.
UnivariatePuiseuxSeriesFunctions2(Coef1,Coef2,var1,var2,cen1,cen2):_
 Exports == Implementation where
  Coef1 : Ring
  Coef2 : Ring
  var1: Symbol
  var2: Symbol
  cen1: Coef1
  cen2: Coef2
  UPS1  ==> UnivariatePuiseuxSeries(Coef1, var1, cen1)
  UPS2  ==> UnivariatePuiseuxSeries(Coef2, var2, cen2)
  ULSP2 ==> UnivariateLaurentSeriesFunctions2(Coef1, Coef2, var1, var2, cen1, cen2)

  Exports ==> with
    map: (Coef1 -> Coef2,UPS1) -> UPS2
      ++ \spad{map(f,g(x))} applies the map f to the coefficients of the
      ++ Puiseux series \spad{g(x)}.

  Implementation ==> add

    map(f,ups) == puiseux(rationalPower ups, map(f, laurentRep ups)$ULSP2)

