here is a preliminary patch, that passes all tests. Sorry about some unrelated
changes (which I won't commit in the same patch, of course)
Martin
Index: combfunc.spad.pamphlet
===================================================================
--- combfunc.spad.pamphlet (revision 393)
+++ combfunc.spad.pamphlet (working copy)
@@ -51,9 +51,6 @@
++ Provides combinatorial functions over an integral domain.
++ Keywords: combinatorial, function, factorial.
++ Examples: )r COMBF INPUT
-
-
-
CombinatorialFunction(R, F): Exports == Implementation where
R: Join(OrderedSet, IntegralDomain)
F: FunctionSpace R
@@ -65,8 +62,6 @@
SMP ==> SparseMultivariatePolynomial(R, K)
Z ==> Integer
- POWER ==> "%power"::Symbol
- OPEXP ==> "exp"::Symbol
SPECIALDIFF ==> "%specialDiff"
SPECIALDISP ==> "%specialDisp"
SPECIALEQUAL ==> "%specialEqual"
@@ -221,7 +216,8 @@
opdsum := operator("%defsum"::Symbol)$CommonOperators
opprod := operator("product"::Symbol)$CommonOperators
opdprod := operator("%defprod"::Symbol)$CommonOperators
- oppow := operator(POWER::Symbol)$CommonOperators
+ oppow := operator("%power"::Symbol)$CommonOperators
+ opexp := operator(operator("exp"::Symbol)$CommonOperators)$F
factorial x == opfact x
binomial(x, y) == opbinom [x, y]
@@ -239,7 +235,7 @@
x ** y ==
-- Do some basic simplifications
- is?(x,POWER) =>
+ is?(x,"%power"::Symbol) =>
args : List F := argument first kernels x
not(#args = 2) => error "Too many arguments to **"
number?(first args) and number?(y) =>
@@ -327,32 +323,14 @@
differentiate(third l, x) * summand l
+ opsum [differentiate(first l, x), second l, third l]
\end{verbatim}
-
-<<package COMBF CombinatorialFunction>>=
- dvprod(l, x) ==
- dm := retract(dummy)@SE
- f := eval(first l, retract(second l)@K, dm::F)
- p := product(f, dm)
-
- opsum [differentiate(first l, x)/first l * p, second l, third l]
-
-
- dvdprod(l, x) ==
- x = retract(y := third l)@SE => 0
- if member?(x, variables(h := third rest rest l)) or
- member?(x, variables(g := third rest l)) then
- error "a product cannot be differentiated with respect to a bound"
- else
- opdsum cons(differentiate(first l, x)/first l, rest l) * opdprod l
-
-@
-The above two operations implement differentiation of products with and without
+The two operations below implement differentiation of products with and without
bounds. Note again, that we cannot even properly define products with bounds
that are not integral.
To differentiate the product, we use Leibniz rule:
$$\frac{d}{dx}\prod_{i=a}^b f(i,x) =
- \sum_{i=a}^b \frac{\frac{d}{dx} f(i,x)}{f(i,x)}\prod_{i=a}^b f(i,x)
+\left(\sum_{i=a}^b \frac{\frac{d}{dx} f(i,x)}{f(i,x)}\right)
+\prod_{i=a}^b f(i,x)
$$
There is one situation where this definition might produce wrong results,
@@ -381,6 +359,25 @@
results. (Issue~\#211)
<<package COMBF CombinatorialFunction>>=
+ dvprod(l, x) ==
+ dm := retract(dummy)@SE
+ f := eval(first l, retract(second l)@K, dm::F)
+ p := product(f, dm)
+
+ opsum [differentiate(first l, x)/first l * p, second l, third l]
+
+
+ dvdprod(l, x) ==
+ x = retract(y := third l)@SE => 0
+ if member?(x, variables(h := third rest rest l)) or
+ member?(x, variables(g := third rest l)) then
+ error "a product cannot be differentiated with respect to a bound"
+ else
+ opdsum cons(differentiate(first l, x)/first l, rest l) * opdprod l
+
+@
+
+<<package COMBF CombinatorialFunction>>=
dprod l ==
prod(summand(l)::O, third(l)::O)
@@ -458,7 +455,7 @@
is?(op, "%defsum"::Symbol) => opdsum
is?(op, "product"::Symbol) => opprod
is?(op, "%defprod"::Symbol) => opdprod
- is?(op, POWER) => oppow
+ is?(op, "%power"::Symbol) => oppow
error "Not a combinatorial operator"
iprod l ==
@@ -517,7 +514,7 @@
*/[eval(first l,k::K,i::F) for i in r1::Z .. r2::Z]
iiipow l ==
- (u := isExpt(x := first l, OPEXP)) case "failed" => kernel(oppow, l)
+ (u := isExpt(x := first l, opexp)) case "failed" => kernel(oppow, l)
rec := u::Record(var: K, exponent: Z)
y := first argument(rec.var)
(r := retractIfCan(y)@Union(Fraction Z, "failed")) case
@@ -540,14 +537,11 @@
zero?(x := first l) =>
zero? second l => error "0 ** 0"
0
--- one? x or zero?(n := second l) => 1
- (x = 1) or zero?(n: F := second l) => 1
--- one? n => x
- (n = 1) => x
- (u := isExpt(x, OPEXP)) case "failed" => kernel(oppow, l)
+ one? x or zero?(n := second l) => 1
+ one? n => x
+ (u := isExpt(x, opexp)) case "failed" => kernel(oppow, l)
rec := u::Record(var: K, exponent: Z)
--- one?(y := first argument(rec.var)) or y = -1 =>
- ((y := first argument(rec.var))=1) or y = -1 =>
+ one?(y := first argument(rec.var)) or y = -1 =>
(operator(rec.var)) (rec.exponent * y * n)
kernel(oppow, l)
Index: transsolve.spad.pamphlet
===================================================================
--- transsolve.spad.pamphlet (revision 393)
+++ transsolve.spad.pamphlet (working copy)
@@ -195,9 +195,11 @@
-- This function was suggested by Manuel Bronstein as a simpler
-- alternative to normalize.
- simplifyingLog(f:RE):RE ==
- (u:=isExpt(f,"exp"::Symbol)) case Record(var:Kernel
RE,exponent:Integer) =>
- rec := u::Record(var:Kernel RE,exponent:Integer)
+ opexp := operator(operator("exp"::Symbol)$CommonOperators)$RE
+ simplifyingLog(f: RE): RE ==
+ (u := isExpt(f, opexp)) case
+ Record(var: Kernel RE, exponent: Integer) =>
+ rec := u::Record(var: Kernel RE, exponent: Integer)
rec.exponent * first argument(rec.var)
log f
Index: op.spad.pamphlet
===================================================================
--- op.spad.pamphlet (revision 393)
+++ op.spad.pamphlet (working copy)
@@ -108,7 +108,7 @@
++ Argument op is modified "in place", i.e. no copy is made.
Implementation ==> add
- -- if narg < 0 then the operator ahs variable arity.
+ -- if narg < 0 then the operator has variable arity.
Rep := Record(opname:Symbol, narg:SingleInteger, props:P)
oper: (Symbol, SingleInteger, P) -> $
Index: fspace.spad.pamphlet
===================================================================
--- fspace.spad.pamphlet (revision 393)
+++ fspace.spad.pamphlet (working copy)
@@ -475,9 +475,6 @@
isExpt:(%,OP) -> Union(Record(var:K,exponent:Z),"failed")
++ isExpt(p,op) returns \spad{[x, n]} if \spad{p = x**n}
++ and \spad{n <> 0} and \spad{x = op(a)}.
- isExpt:(%,SY) -> Union(Record(var:K,exponent:Z),"failed")
- ++ isExpt(p,f) returns \spad{[x, n]} if \spad{p = x**n}
- ++ and \spad{n <> 0} and \spad{x = f(a)}.
isPower : % -> Union(Record(val:%,exponent:Z),"failed")
++ isPower(p) returns \spad{[x, n]} if \spad{p = x**n}
++ and \spad{n <> 0}.
@@ -831,16 +828,9 @@
isExpt(x:%, op:OP) ==
(u := isExpt x) case "failed" => "failed"
- v := (u::Record(var:K, exponent:Z)).var
- is?(v, op) and #argument(v) = 1 => u
+ is?((u::Record(var:K, exponent:Z)).var, op) => u
"failed"
- isExpt(x:%, sy:SY) ==
- (u := isExpt x) case "failed" => "failed"
- v := (u::Record(var:K, exponent:Z)).var
- is?(v, sy) and #argument(v) = 1 => u
- "failed"
-
if R has RetractableTo Z then
smpIsMult p ==
-- (u := mainVariable p) case K and one?
degree(q:=univariate(p,u::K))
Index: algfunc.spad.pamphlet
===================================================================
--- algfunc.spad.pamphlet (revision 393)
+++ algfunc.spad.pamphlet (working copy)
@@ -440,11 +440,16 @@
x ** qr.quotient * inroot([x, (denom q)::F]) ** qr.remainder
<<hackroot(x, n)>>
+@
+The following defines the evaluation function for \spad{oproot}, which has the
+name \spad{nthRoot}. It's first argument is the radicand, the second the \lq
+root-exponent\rq, i.e. $\sqrt[n]{x}$ is represented as \spad{nthRoot [x, n]}.
+
+<<package AF AlgebraicFunction>>=
inroot l ==
zero?(n := retract(second l)@Z) => error "root: exponent = 0"
--- one?(x := first l) or one? n => x
- ((x := first l) = 1) or (n = 1) => x
+ one?(x := first l) or one? n => x
(r := retractIfCan(x)@Union(R,"failed")) case R => iroot(r::R,n)
(u := isExpt(x)) case Record(var:K, exponent:Z) =>
pr := u::Record(var:K, exponent:Z)
--~--~---------~--~----~------------~-------~--~----~
You received this message because you are subscribed to the Google Groups
"FriCAS - computer algebra system" group.
To post to this group, send email to [email protected]
To unsubscribe from this group, send email to [EMAIL PROTECTED]
For more options, visit this group at
http://groups.google.com/group/fricas-devel?hl=en
-~----------~----~----~----~------~----~------~--~---
