login  home  contents  what's new  discussion  bug reports     help  links  subscribe  changes  refresh  edit

Edit detail for CelineMethod revision 1 of 3

1 2 3
Editor: pagani
Time: 2026/05/31 15:04:15 GMT+0
Note:

changed:
-
\begin{spad}
)abbrev package CELINE CelineMethod
CelineMethod: Exports == Implementation where
  -- Ref: https://www2.math.upenn.edu/~wilf/AeqB.pdf
  -- R:Join(Comparable,IntegralDomain,RTR,LER) 
  R   ==> Integer  
  X   ==> Expression R
  FPR ==> Fraction Polynomial R
  
  RTR ==> RetractableTo R
  LER ==> LinearlyExplicitOver R
  
  POLR ==> Polynomial R
  NNI  ==> NonNegativeInteger
  PI   ==> PositiveInteger
  BOP  ==> BasicOperator
  
  EFSP ==> ElementaryFunctionStructurePackage(Integer,X)
  UPR  ==> Union(Polynomial R,"failed")
  SMP  ==> SparseMultivariatePolynomial(Integer,Kernel(X))
  SMPF ==> SparseMultivariatePolynomial(R,Kernel(X))
  FSMP ==> Factored SMP



  Exports == with
    coeffMatrix : (PI,PI) -> Matrix X
      ++ coeffMatrix(I,J) constructs a (I+1)x(J+1) matrix of new
      ++ symbols.
    funcMatrix  : ((X,X)->X,PI,PI) -> Matrix X
      ++ funcMatrix(f,I,J) constructs a (J+1)x(I+1) matrix with
      ++ the entries M(i,j):=f(n-j,k-i)
    sumAndNormalize : (Matrix X,Matrix X) -> X
      ++ 
    retractNumeratorToPolyInt : X -> UPR
      ++
    findCoeffs : (POLR,Matrix X, Symbol) -> List List Equation FPR
      ++
    celine :((X,X)->X,PI,PI) -> Union(Matrix X,"failed")
      ++ celine(f,I,J) tries to find a recurrence for 
      ++ F(n) = summation(f(n,k),k). The returned value is either
      ++ "failed" or the coefficient matrix for the recurrence.
      ++ If successfull, then it holds that the trace of 
      ++ celine(f,I,J)*funcMatrix(f,I,J) is zero (normalize). 
    firstOrderRecurrence : ((X,X)->X,X->Matrix X) -> X 
      ++
    recRel : ((X,X)->X,PI,PI) -> Union(X,"failed")  
      ++
    recEq : ((X,X)->X,PI,PI) -> Union(Equation X,"failed")          
      ++ 


  Implementation == add

    coeffMatrix(I:PI,J:PI):Matrix X ==
      matrix [[new()$Symbol::X for j in 0..J] for i in 0..I]

    funcMatrix(F:(X,X)->X,I:PI,J:PI):Matrix X ==
      n:X:='n::X
      k:X:='k::X
      matrix [[F(n-j::X,k-i::X)/F(n,k)  for i in 0..I] for j in 0..J]

    sumAndNormalize(CM:Matrix X, FM:Matrix X):X ==
      n:=nrows(CM)
      t:X:=trace(squareMatrix(CM*FM)$SquareMatrix(n,X))
      normalize(t)$ElementaryFunctionStructurePackage(Integer,X)

    retractNumeratorToPolyInt(san:X):UPR ==
      --p:Polynomial R:=retract(numerator san)
      p:=retractIfCan(numerator san)@UPR
      p case "failed" => "failed"
      return p
 

    findCoeffs(p:POLR,m:Matrix X,k:Symbol):List List Equation FPR ==
      d:=degree(p,k)
      eqs:List Equation FPR:=[coefficient(p,k,l)::FPR=0 for l in 0..d]
      v:=members m
      x:=variables v
      sol:=solve(eqs,x)$SystemSolvePackage(R) --$TransSolvePackage(Integer)
      --first sol -- if more?
      
      
    convToFPR(M:Matrix X):Matrix FPR ==
      m:=copy(M)
      nr:=nrows(m)
      nc:=ncols(m)
      r:=zero(nr,nc)$Matrix(FPR)
      for i in 1..nr repeat
        for j in 1..nc repeat
          r(i,j):=retract m(i,j)
      return r

    convertToX(x:List Equation FPR):List Equation X ==
      [lhs(s)::X=rhs(s)::X for s in x]

    celine(F:(X,X)->X,I:PI,J:PI):Union(Matrix X,"failed") ==
      cm:Matrix X:=coeffMatrix(I,J)
      fm:=funcMatrix(F,I,J)
      san:=sumAndNormalize(cm,fm)
      p:=retractNumeratorToPolyInt(san)
      p case "failed" => "failed"
      --
      fseq:=first findCoeffs(p,cm,'k)
      seq:=convertToX(fseq)
      ev:=(x:X):X+->eval(x,seq)
      map(ev,cm)


    recRel(F:(X,X)->X,I:PI,J:PI):Union(X,"failed") ==
      C:=celine(F,I,J)
      C case "failed" => "failed"
      zero? C => "failed"
      S:=operator 'S
      SM:=funcMatrix((n,k)+->S(n),I,J)
      if C case Matrix(X) then
        M:=C*SM      
        r:=numerator reduce(_+,[M(i,i) for i in 1..I+1])
      else
        "failed"

    -- solve([r=0],[S(n)])
    -- F:=(n,k)+->binomial(n,k)*binomial(2*k,k)*(-2)^(n-k)
    
    recEq(F:(X,X)->X,I:PI,J:PI):Union(Equation X,"failed") ==
      r:Union(X,"failed"):=recRel(F,I,J)
      r case "failed" => "failed"
      S:=operator 'S
      n:='n::X
      s:=solve([r=0],[S(n)::X])$TransSolvePackage(Integer)
      empty? s => "failed"
      first first s

    firstOrderRecurrence(F:(X,X)->X,A:X->Matrix X):X ==
      n:X:='n::X
      k:X:='k::X
      --fm:=funcMatrix(F,1,1)
      a00:=(n:X):X+->elt(A(n),1,1)
      a01:=(n:X):X+->elt(A(n),1,2)
      a10:=(n:X):X+->elt(A(n),2,1)
      a11:=(n:X):X+->elt(A(n),2,2)
      anum:=(n:X):X+->a11(n)+a01(n)
      aden:=(n:X):X+->a00(n)+a10(n)
      alpha:=(n:X):X+->-anum(n)/aden(n)
      beta:=(n:X):X+->a10(n)/aden(n)*F(n+1,n+1)
      sb1:SegmentBinding X:=equation('k,segment(0,n)$Segment(X))
      sb2:SegmentBinding X:=equation('k,segment(1,n)$Segment(X))
      Q:=(n:X):X+->product(alpha(k),sb1)
      if F(0,0) ~= 0 then
        r:X:=Q(n)*(F(0,0)/alpha(0)+ summation(beta(k)/Q(k),sb2))
      else
        r:X:=Q(n)*summation(beta(k)/Q(k),sb2)
      return r
      

-- https://en.wikipedia.org/wiki/Recurrence_relation

\end{spad}

\begin{axiom}

)version

X ==> EXPR INT


F1:=(n,k)+->binomial(n,k)
r1:=recEq(F1,1,1)

-- Example 4.3.1. 
F2:=(n,k)+->binomial(n,k)^2
r2:=recEq(F2,2,2)

-- Example 4.3.3. 
F3:=(n,k)+->binomial(n,k)*binomial(2*k,k)*(-2)^(n-k)
r3:=recEq(F3,2,2)

-- Example 4.3.2. 
F4:=(n,k)+->(-1)^k*binomial(n,k)*x^k/factorial(k)
r4:=recEq(F4,2,2)  -- (1,2) already works!

-- Example 4.1.1.
F5:=(n,k)+->k*binomial(n,k)
r5:=recEq(F5,1,1)


F6:=(n,k)+->k^2*binomial(n,k)
r6:=recEq(F6,1,1) -- fails, but (1,2) and (2,1) are working .... pref: 1,2

r6_12:=recEq(F6,1,2)
r6_21:=recEq(F6,2,1)

r6_21 - r612

\end{axiom}


fricas
(1) -> <spad>
fricas
)abbrev package CELINE CelineMethod
CelineMethod: Exports == Implementation where
  -- Ref: https://www2.math.upenn.edu/~wilf/AeqB.pdf
  -- R:Join(Comparable,IntegralDomain,RTR,LER) 
  R   ==> Integer  
  X   ==> Expression R
  FPR ==> Fraction Polynomial R
RTR ==> RetractableTo R LER ==> LinearlyExplicitOver R
POLR ==> Polynomial R NNI ==> NonNegativeInteger PI ==> PositiveInteger BOP ==> BasicOperator
EFSP ==> ElementaryFunctionStructurePackage(Integer,X) UPR ==> Union(Polynomial R,"failed") SMP ==> SparseMultivariatePolynomial(Integer,Kernel(X)) SMPF ==> SparseMultivariatePolynomial(R,Kernel(X)) FSMP ==> Factored SMP
Exports == with coeffMatrix : (PI,PI) -> Matrix X ++ coeffMatrix(I,J) constructs a (I+1)x(J+1) matrix of new ++ symbols. funcMatrix : ((X,X)->X,PI,PI) -> Matrix X ++ funcMatrix(f,I,J) constructs a (J+1)x(I+1) matrix with ++ the entries M(i,j):=f(n-j,k-i) sumAndNormalize : (Matrix X,Matrix X) -> X ++ retractNumeratorToPolyInt : X -> UPR ++ findCoeffs : (POLR,Matrix X, Symbol) -> List List Equation FPR ++ celine :((X,X)->X,PI,PI) -> Union(Matrix X,"failed") ++ celine(f,I,J) tries to find a recurrence for ++ F(n) = summation(f(n,k),k). The returned value is either ++ "failed" or the coefficient matrix for the recurrence. ++ If successfull, then it holds that the trace of ++ celine(f,I,J)*funcMatrix(f,I,J) is zero (normalize). firstOrderRecurrence : ((X,X)->X,X->Matrix X) -> X ++ recRel : ((X,X)->X,PI,PI) -> Union(X,"failed") ++ recEq : ((X,X)->X,PI,PI) -> Union(Equation X,"failed") ++
Implementation == add
coeffMatrix(I:PI,J:PI):Matrix X == matrix [[new()$Symbol::X for j in 0..J] for i in 0..I]
funcMatrix(F:(X,X)->X,I:PI,J:PI):Matrix X == n:X:='n::X k:X:='k::X matrix [[F(n-j::X,k-i::X)/F(n,k) for i in 0..I] for j in 0..J]
sumAndNormalize(CM:Matrix X, FM:Matrix X):X == n:=nrows(CM) t:X:=trace(squareMatrix(CM*FM)$SquareMatrix(n,X)) normalize(t)$ElementaryFunctionStructurePackage(Integer,X)
retractNumeratorToPolyInt(san:X):UPR == --p:Polynomial R:=retract(numerator san) p:=retractIfCan(numerator san)@UPR p case "failed" => "failed" return p
findCoeffs(p:POLR,m:Matrix X,k:Symbol):List List Equation FPR == d:=degree(p,k) eqs:List Equation FPR:=[coefficient(p,k,l)::FPR=0 for l in 0..d] v:=members m x:=variables v sol:=solve(eqs,x)$SystemSolvePackage(R) --$TransSolvePackage(Integer) --first sol -- if more?
convToFPR(M:Matrix X):Matrix FPR == m:=copy(M) nr:=nrows(m) nc:=ncols(m) r:=zero(nr,nc)$Matrix(FPR) for i in 1..nr repeat for j in 1..nc repeat r(i,j):=retract m(i,j) return r
convertToX(x:List Equation FPR):List Equation X == [lhs(s)::X=rhs(s)::X for s in x]
celine(F:(X,X)->X,I:PI,J:PI):Union(Matrix X,"failed") == cm:Matrix X:=coeffMatrix(I,J) fm:=funcMatrix(F,I,J) san:=sumAndNormalize(cm,fm) p:=retractNumeratorToPolyInt(san) p case "failed" => "failed" -- fseq:=first findCoeffs(p,cm,'k) seq:=convertToX(fseq) ev:=(x:X):X+->eval(x,seq) map(ev,cm)
recRel(F:(X,X)->X,I:PI,J:PI):Union(X,"failed") == C:=celine(F,I,J) C case "failed" => "failed" zero? C => "failed" S:=operator 'S SM:=funcMatrix((n,k)+->S(n),I,J) if C case Matrix(X) then M:=C*SM r:=numerator reduce(_+,[M(i,i) for i in 1..I+1]) else "failed"
-- solve([r=0],[S(n)]) -- F:=(n,k)+->binomial(n,k)*binomial(2*k,k)*(-2)^(n-k)
recEq(F:(X,X)->X,I:PI,J:PI):Union(Equation X,"failed") == r:Union(X,"failed"):=recRel(F,I,J) r case "failed" => "failed" S:=operator 'S n:='n::X s:=solve([r=0],[S(n)::X])$TransSolvePackage(Integer) empty? s => "failed" first first s
firstOrderRecurrence(F:(X,X)->X,A:X->Matrix X):X == n:X:='n::X k:X:='k::X --fm:=funcMatrix(F,1,1) a00:=(n:X):X+->elt(A(n),1,1) a01:=(n:X):X+->elt(A(n),1,2) a10:=(n:X):X+->elt(A(n),2,1) a11:=(n:X):X+->elt(A(n),2,2) anum:=(n:X):X+->a11(n)+a01(n) aden:=(n:X):X+->a00(n)+a10(n) alpha:=(n:X):X+->-anum(n)/aden(n) beta:=(n:X):X+->a10(n)/aden(n)*F(n+1,n+1) sb1:SegmentBinding X:=equation('k,segment(0,n)$Segment(X)) sb2:SegmentBinding X:=equation('k,segment(1,n)$Segment(X)) Q:=(n:X):X+->product(alpha(k),sb1) if F(0,0) ~= 0 then r:X:=Q(n)*(F(0,0)/alpha(0)+ summation(beta(k)/Q(k),sb2)) else r:X:=Q(n)*summation(beta(k)/Q(k),sb2) return r
-- https://en.wikipedia.org/wiki/Recurrence_relation</spad>
fricas
Compiling FriCAS source code from file 
      /var/lib/zope2.10/instance/axiom-wiki/var/LatexWiki/5910754882400853191-25px001.spad
      using old system compiler.
   CELINE abbreviates package CelineMethod 
------------------------------------------------------------------------
   initializing NRLIB CELINE for CelineMethod 
   compiling into NRLIB CELINE 
   compiling exported coeffMatrix : (PositiveInteger,PositiveInteger) -> Matrix Expression Integer
Time: 0.03 SEC.
compiling exported funcMatrix : ((Expression Integer,Expression Integer) -> Expression Integer,PositiveInteger,PositiveInteger) -> Matrix Expression Integer Time: 0.01 SEC.
compiling exported sumAndNormalize : (Matrix Expression Integer,Matrix Expression Integer) -> Expression Integer Time: 0.02 SEC.
compiling exported retractNumeratorToPolyInt : Expression Integer -> Union(Polynomial Integer,failed) Time: 0 SEC.
compiling exported findCoeffs : (Polynomial Integer,Matrix Expression Integer,Symbol) -> List List Equation Fraction Polynomial Integer Time: 0.01 SEC.
compiling local convToFPR : Matrix Expression Integer -> Matrix Fraction Polynomial Integer Time: 0 SEC.
compiling local convertToX : List Equation Fraction Polynomial Integer -> List Equation Expression Integer Time: 0 SEC.
compiling exported celine : ((Expression Integer,Expression Integer) -> Expression Integer,PositiveInteger,PositiveInteger) -> Union(Matrix Expression Integer,failed) Time: 0 SEC.
compiling exported recRel : ((Expression Integer,Expression Integer) -> Expression Integer,PositiveInteger,PositiveInteger) -> Union(Expression Integer,failed) Time: 0.04 SEC.
compiling exported recEq : ((Expression Integer,Expression Integer) -> Expression Integer,PositiveInteger,PositiveInteger) -> Union(Equation Expression Integer,failed) Time: 0.06 SEC.
compiling exported firstOrderRecurrence : ((Expression Integer,Expression Integer) -> Expression Integer,Expression Integer -> Matrix Expression Integer) -> Expression Integer Time: 0.20 SEC.
(time taken in buildFunctor: 0) Time: 0 SEC.
Warnings: [1] sumAndNormalize: not known that (AlgebraicallyClosedField) is of mode (CATEGORY domain (IF (has (Integer) (IntegralDomain)) (PROGN (ATTRIBUTE (AlgebraicallyClosedFunctionSpace (Integer))) (ATTRIBUTE (TranscendentalFunctionCategory)) (ATTRIBUTE (CombinatorialOpsCategory)) (ATTRIBUTE (LiouvillianFunctionCategory)) (ATTRIBUTE (SpecialFunctionCategory)) (SIGNATURE reduce (% %)) (SIGNATURE number? ((Boolean) %)) (IF (has (Integer) (PolynomialFactorizationExplicit)) (ATTRIBUTE (PolynomialFactorizationExplicit)) noBranch) (SIGNATURE setSimplifyDenomsFlag ((Boolean) (Boolean))) (SIGNATURE getSimplifyDenomsFlag ((Boolean)))) noBranch)) [2] sumAndNormalize: not known that (TranscendentalFunctionCategory) is of mode (CATEGORY domain (IF (has (Integer) (IntegralDomain)) (PROGN (ATTRIBUTE (AlgebraicallyClosedFunctionSpace (Integer))) (ATTRIBUTE (TranscendentalFunctionCategory)) (ATTRIBUTE (CombinatorialOpsCategory)) (ATTRIBUTE (LiouvillianFunctionCategory)) (ATTRIBUTE (SpecialFunctionCategory)) (SIGNATURE reduce (% %)) (SIGNATURE number? ((Boolean) %)) (IF (has (Integer) (PolynomialFactorizationExplicit)) (ATTRIBUTE (PolynomialFactorizationExplicit)) noBranch) (SIGNATURE setSimplifyDenomsFlag ((Boolean) (Boolean))) (SIGNATURE getSimplifyDenomsFlag ((Boolean)))) noBranch))
Cumulative Statistics for Constructor CelineMethod Time: 0.39 seconds
finalizing NRLIB CELINE Processing CelineMethod for Browser database: --->-->CelineMethod(constructor): Not documented!!!! --------(coeffMatrix ((Matrix (Expression (Integer))) (PositiveInteger) (PositiveInteger)))--------- --------(funcMatrix ((Matrix (Expression (Integer))) (Mapping (Expression (Integer)) (Expression (Integer)) (Expression (Integer))) (PositiveInteger) (PositiveInteger)))--------- --------(sumAndNormalize ((Expression (Integer)) (Matrix (Expression (Integer))) (Matrix (Expression (Integer)))))--------- --------(retractNumeratorToPolyInt ((Union (Polynomial (Integer)) failed) (Expression (Integer))))--------- --------(findCoeffs ((List (List (Equation (Fraction (Polynomial (Integer)))))) (Polynomial (Integer)) (Matrix (Expression (Integer))) (Symbol)))--------- --------(celine ((Union (Matrix (Expression (Integer))) failed) (Mapping (Expression (Integer)) (Expression (Integer)) (Expression (Integer))) (PositiveInteger) (PositiveInteger)))--------- --------(firstOrderRecurrence ((Expression (Integer)) (Mapping (Expression (Integer)) (Expression (Integer)) (Expression (Integer))) (Mapping (Matrix (Expression (Integer))) (Expression (Integer)))))--------- --------(recRel ((Union (Expression (Integer)) failed) (Mapping (Expression (Integer)) (Expression (Integer)) (Expression (Integer))) (PositiveInteger) (PositiveInteger)))--------- --------(recEq ((Union (Equation (Expression (Integer))) failed) (Mapping (Expression (Integer)) (Expression (Integer)) (Expression (Integer))) (PositiveInteger) (PositiveInteger)))--------- --->-->CelineMethod(): Missing Description ; compiling file "/var/aw/var/LatexWiki/CELINE.NRLIB/CELINE.lsp" (written 31 MAY 2026 03:04:16 PM):
; wrote /var/aw/var/LatexWiki/CELINE.NRLIB/CELINE.fasl ; compilation finished in 0:00:00.068 ------------------------------------------------------------------------ CelineMethod is now explicitly exposed in frame initial CelineMethod will be automatically loaded when needed from /var/aw/var/LatexWiki/CELINE.NRLIB/CELINE

fricas
)version
"FriCAS 1.3.12 compiled at Sat 7 Jun 23:54:49 CEST 2025"
X ==> EXPR INT
Type: Void
fricas
F1:=(n,k)+->binomial(n,k)

\label{eq1}{\left(n , \: k \right)}\mapsto{binomial \left({n , \: k}\right)}(1)
Type: AnonymousFunction?
fricas
r1:=recEq(F1,1,1)

\label{eq2}{S \left({n}\right)}={2 \ {S \left({n - 1}\right)}}(2)
Type: Union(Equation(Expression(Integer)),...)
fricas
-- Example 4.3.1. 
F2:=(n,k)+->binomial(n,k)^2

\label{eq3}{\left(n , \: k \right)}\mapsto{{binomial \left({n , \: k}\right)}^{2}}(3)
Type: AnonymousFunction?
fricas
r2:=recEq(F2,2,2)

\label{eq4}{S \left({n}\right)}={\frac{{\left({4 \  n}- 2 \right)}\ {S \left({n - 1}\right)}}{n}}(4)
Type: Union(Equation(Expression(Integer)),...)
fricas
-- Example 4.3.3. 
F3:=(n,k)+->binomial(n,k)*binomial(2*k,k)*(-2)^(n-k)

\label{eq5}{\left(n , \: k \right)}\mapsto{{binomial \left({n , \: k}\right)}\ {binomial \left({{2 \  k}, \: k}\right)}\ {{\left(- 2 \right)}^{n - k}}}(5)
Type: AnonymousFunction?
fricas
r3:=recEq(F3,2,2)

\label{eq6}{S \left({n}\right)}={\frac{{\left({4 \  n}- 4 \right)}\ {S \left({n - 2}\right)}}{n}}(6)
Type: Union(Equation(Expression(Integer)),...)
fricas
-- Example 4.3.2. 
F4:=(n,k)+->(-1)^k*binomial(n,k)*x^k/factorial(k)

\label{eq7}{\left(n , \: k \right)}\mapsto{\frac{{{\left(- 1 \right)}^{k}}\ {binomial \left({n , \: k}\right)}\ {{x}^{k}}}{factorial \left({k}\right)}}(7)
Type: AnonymousFunction?
fricas
r4:=recEq(F4,2,2)  -- (1,2) already works!

\label{eq8}{S \left({n}\right)}={\frac{{{\left(- x +{2 \  n}- 1 \right)}\ {S \left({n - 1}\right)}}+{{\left(- n + 1 \right)}\ {S \left({n - 2}\right)}}}{n}}(8)
Type: Union(Equation(Expression(Integer)),...)
fricas
-- Example 4.1.1.
F5:=(n,k)+->k*binomial(n,k)

\label{eq9}{\left(n , \: k \right)}\mapsto{k \ {binomial \left({n , \: k}\right)}}(9)
Type: AnonymousFunction?
fricas
r5:=recEq(F5,1,1)

\label{eq10}{S \left({n}\right)}={\frac{2 \  n \ {S \left({n - 1}\right)}}{n - 1}}(10)
Type: Union(Equation(Expression(Integer)),...)
fricas
F6:=(n,k)+->k^2*binomial(n,k)

\label{eq11}{\left(n , \: k \right)}\mapsto{{{k}^{2}}\ {binomial \left({n , \: k}\right)}}(11)
Type: AnonymousFunction?
fricas
r6:=recEq(F6,1,1) -- fails, but (1,2) and (2,1) are working .... pref: 1,2

\label{eq12}\verb#"failed"#(12)
Type: Union("failed",...)
fricas
r6_12:=recEq(F6,1,2)

\label{eq13}{S \left({n}\right)}={\frac{{3 \  n \ {S \left({n - 1}\right)}}-{2 \  n \ {S \left({n - 2}\right)}}}{n - 1}}(13)
Type: Union(Equation(Expression(Integer)),...)
fricas
r6_21:=recEq(F6,2,1)

\label{eq14}{S \left({n}\right)}={\frac{{\left({2 \  n}+ 2 \right)}\ {S \left({n - 1}\right)}}{n - 1}}(14)
Type: Union(Equation(Expression(Integer)),...)
fricas
r6_21 - r612

\label{eq15}{{S \left({n}\right)}- r 612}={\frac{{{\left({2 \  n}+ 2 \right)}\ {S \left({n - 1}\right)}}+{{\left(- n + 1 \right)}\  r 612}}{n - 1}}(15)
Type: Equation(Expression(Integer))